Statistical Analysis
Jan Taubenheim
7/30/2020
knitr::opts_chunk$set(message = FALSE, warning = FALSE)
require(data.table)
require(DESeq2)
if(!("qiime2R" %in% installed.packages()[,"Package"])){
devtools::install_github("jbisanz/qiime2R")
}
require(qiime2R)
require(ggplot2)
require(vegan)
require(phyloseq)
require(umap)
require(MASS)
require(cowplot)
require(geosphere)
require(lme4)
require(nlme)
require(car)
require(jtools)Statistical analysis
The main questions are:
* do microbiomes of Hydra from different locations differ
* do microbiomes differ between species
* do microbiomes differ between different sexual states The two common ways to test that is using an Anosim (a non parametric test on beta diversity) and an Adonis test (comparable to MANOVA on beta diversity). Common practice is to resample the libraries to the minimum number of reads per sample with the problem of loosing most of the data. However, the recommendation for statistical analysis based on the phyloseq authors are to normalize the reads with the DESeq2 package (which is also used to determine differential ASV between conditions) and the implementation of the variance stabilization transformation (vst). This requires, that at least one ASV is common to all samples, which is not the case for the diverse sample I got. To circumvent this problem, I will add a single read to every ASV in every sample. Another common and quite well working transformation is dividing by the library size. I will use the raw, the vst and the fraction to do the statistics.
otu <- read.csv(snakemake@input[["AbFilt_FT"]], row.names = 1) #"../data/dada2_AbFilt_FT.csv"
tax <- read.csv(snakemake@input[["AbFilt_Tax"]], row.names=1)#"../data/dada2_AbFilt_Tax.csv"
meta <- read.csv(snakemake@input[["AbFilt_meta"]], row.names=1)#"../data/dada2_AbFilt_meta.csv"
tree <- read_qza(snakemake@input[["AbFilt_fasta_tree"]])$data #"../data/dada2_AbFilt_fasta_tree.qza"
#otu <- read.csv("../data/dada2_AbFilt_FT.csv", row.names = 1)
#tax <- read.csv("../data/dada2_AbFilt_Tax.csv", row.names =1)
#meta <- read.csv("../data/dada2_AbFilt_meta.csv", row.names=1)
#tree <- read_qza("../data/dada2_AbFilt_fasta_tree.qza")$data
# remap the reproductive mode to the SEX, ASEX and NR
map_vect <- c("SEX","NR","SEX","ASEX","SEX","SEX","SEX","SEX")
names(map_vect) <- unique(meta$ReproductiveMode)
meta[["ReproductiveModeSimple"]] <- map_vect[meta[["ReproductiveMode"]]]
# artificially add a single count to all ESVs
otu2 <- otu+1
phy <- phyloseq(otu_table(otu, taxa_are_rows = TRUE),
tax_table(as.matrix(tax)),
sample_data(meta),
phy_tree(tree))
otu_vst <- varianceStabilizingTransformation(as.matrix(otu2), fitType = "local")
# shift the vst correction to positive values
if(min(otu_vst) <0){
otu_vst <- otu_vst + min(otu_vst)*-1
}
# create phyloseq objects
phy_vst <- phyloseq(otu_table(otu_vst,
taxa_are_rows = TRUE),
tax_table(as.matrix(tax)),
sample_data(meta),
phy_tree(tree))
phy_frac <- phyloseq(otu_table(otu/colSums(otu), taxa_are_rows = TRUE),
tax_table(as.matrix(tax)),
sample_data(meta),
phy_tree(tree))
trnses <- list(raw = phy, vst = phy_vst, frac = phy_frac)So lets look on the raw data here and whether we already see a separation there (on euclidean distance).
countPlots <- list()
plotDF <- data.frame()
for(n in names(trnses)){
# calculate ordination
pca <- prcomp(t(otu_table(trnses[[n]])))
nmds <- isoMDS(dist(t(otu_table(trnses[[n]]))))
ump <- umap(t(otu_table(trnses[[n]])))
# create plots
nPCA <- paste(n,"PCA",sep = "_")
countPlots[[nPCA]] <- ggplot(as.data.frame(pca$x),
aes(x=PC1, y=PC2, col = meta$PopID)) +
geom_point() +
ggtitle(nPCA) +
theme(legend.position = "none")
nNMDS <- paste(n,"NMDS", sep = "_")
countPlots[[nNMDS]] <- ggplot(as.data.frame(nmds$points),
aes(x=V1, y=V2, col = meta$PopID)) +
geom_point() +
ggtitle(nNMDS) +
theme(legend.position = "none")
nUMAP <- paste(n,"UMAP",sep ="_")
countPlots[[nUMAP]] <- ggplot(as.data.frame(ump$layout),
aes(x=V1, y=V2, col = meta$PopID))+
geom_point() +
ggtitle(nUMAP)
# combine the plotting data in one data.frame
plotDF <- rbind(plotDF,
cbind(V1 = pca$x[,1], V2 =pca$x[,2], ordination = "PCA", data = n, sample = rownames(pca$x)),
cbind(as.data.frame(nmds$points), ordination = "NMDS", data = n, sample = rownames(nmds$points)),
cbind(as.data.frame(ump$layout), ordination = "UMAP", data =n, sample = rownames(ump$layout)),
deparse.level = 2
)
}## initial value 39.315036
## final value 39.315035
## converged
## initial value 37.542053
## final value 37.523047
## converged
## initial value 54.716981
## iter 5 value 13.867358
## iter 10 value 12.441360
## iter 15 value 12.212251
## iter 20 value 12.078578
## iter 20 value 12.071296
## iter 20 value 12.071248
## final value 12.071248
## converged# plot the different normalizations next to each other
plotDF<- merge(plotDF, meta, by.x="sample", by.y=0)
plotDF[,"V1"] <- as.numeric(plotDF[,"V1"])
plotDF[,"V2"] <- as.numeric(plotDF[,"V2"])
normPlots <- list()
for(n in names(trnses)){
normPlots[[n]] <- ggplot(plotDF[plotDF[["data"]] == n,], aes(x = V1, y=V2, col =PopID)) +
facet_wrap(.~ordination, scales = "free") +
geom_point()+
ggtitle(n) +
theme(legend.position = "None")
}
cowplot::plot_grid(plotlist = normPlots, ncol =1)
To calculate differences between treatments, one needs to calculate \(\beta\)-diversity between samples. There are several distance measures, which can be calculated, the most common though are: Bray-Curtis, Jaccard, Unifrac and Weighted Unifrac. I will calculate all. To evaluate which might be the best to do statistics on, and which input data to use, I will plot ordination of the distance data. This is done by NMDS, PCoA, and UMAP.
## [1] "euclidean_raw"
## [1] "bray_raw"
## [1] "jaccard_raw"
## [1] "unifrac_raw"
## [1] "wunifrac_raw"
## [1] "euclidean_vst"
## [1] "bray_vst"
## [1] "jaccard_vst"
## [1] "unifrac_vst"
## [1] "wunifrac_vst"
## [1] "euclidean_frac"
## [1] "bray_frac"
## [1] "jaccard_frac"
## [1] "unifrac_frac"
## [1] "wunifrac_frac"## [1] "euclidean_raw_NMDS"
## Run 0 stress 0.08369916
## Run 1 stress 0.08493086
## Run 2 stress 0.08730913
## Run 3 stress 0.08644224
## Run 4 stress 0.08512591
## Run 5 stress 0.08413704
## ... Procrustes: rmse 0.03578869 max resid 0.3719535
## Run 6 stress 0.08379216
## ... Procrustes: rmse 0.0224323 max resid 0.1555071
## Run 7 stress 0.08985744
## Run 8 stress 0.08678246
## Run 9 stress 0.08613573
## Run 10 stress 0.08426065
## Run 11 stress 0.0838994
## ... Procrustes: rmse 0.02361745 max resid 0.2184259
## Run 12 stress 0.08653453
## Run 13 stress 0.08513959
## Run 14 stress 0.0867496
## Run 15 stress 0.08454579
## Run 16 stress 0.08578121
## Run 17 stress 0.08490396
## Run 18 stress 0.08705293
## Run 19 stress 0.08369116
## ... New best solution
## ... Procrustes: rmse 0.02715217 max resid 0.2166256
## Run 20 stress 0.08508169
## *** No convergence -- monoMDS stopping criteria:
## 18: no. of iterations >= maxit
## 2: stress ratio > sratmax
## [1] "euclidean_raw_PCoA"
## [1] "euclidean_raw_UMAP"
## [1] "bray_raw_NMDS"
## Run 0 stress 0.2773178
## Run 1 stress 0.2811224
## Run 2 stress 0.2814037
## Run 3 stress 0.2797893
## Run 4 stress 0.2781847
## Run 5 stress 0.278967
## Run 6 stress 0.2861727
## Run 7 stress 0.2798618
## Run 8 stress 0.2770567
## ... New best solution
## ... Procrustes: rmse 0.02396036 max resid 0.1512107
## Run 9 stress 0.2769995
## ... New best solution
## ... Procrustes: rmse 0.02232634 max resid 0.1157627
## Run 10 stress 0.2789943
## Run 11 stress 0.2845254
## Run 12 stress 0.2773256
## ... Procrustes: rmse 0.0263711 max resid 0.1499273
## Run 13 stress 0.3120813
## Run 14 stress 0.2812445
## Run 15 stress 0.2772958
## ... Procrustes: rmse 0.02445731 max resid 0.1212653
## Run 16 stress 0.2807581
## Run 17 stress 0.277671
## Run 18 stress 0.282163
## Run 19 stress 0.2783355
## Run 20 stress 0.2779285
## *** No convergence -- monoMDS stopping criteria:
## 5: no. of iterations >= maxit
## 15: stress ratio > sratmax
## [1] "bray_raw_PCoA"
## [1] "bray_raw_UMAP"
## [1] "jaccard_raw_NMDS"
## Run 0 stress 0.2773784
## Run 1 stress 0.2789441
## Run 2 stress 0.2787378
## Run 3 stress 0.3166571
## Run 4 stress 0.2806936
## Run 5 stress 0.280006
## Run 6 stress 0.2770098
## ... New best solution
## ... Procrustes: rmse 0.01917113 max resid 0.1242275
## Run 7 stress 0.2842891
## Run 8 stress 0.280299
## Run 9 stress 0.2880282
## Run 10 stress 0.2812404
## Run 11 stress 0.2781231
## Run 12 stress 0.2786286
## Run 13 stress 0.2802902
## Run 14 stress 0.2769373
## ... New best solution
## ... Procrustes: rmse 0.02350488 max resid 0.155115
## Run 15 stress 0.2806287
## Run 16 stress 0.2821839
## Run 17 stress 0.2811259
## Run 18 stress 0.2792081
## Run 19 stress 0.278262
## Run 20 stress 0.2759648
## ... New best solution
## ... Procrustes: rmse 0.01679744 max resid 0.1232613
## *** No convergence -- monoMDS stopping criteria:
## 8: no. of iterations >= maxit
## 11: stress ratio > sratmax
## 1: scale factor of the gradient < sfgrmin
## [1] "jaccard_raw_PCoA"
## [1] "jaccard_raw_UMAP"
## [1] "unifrac_raw_NMDS"
## Run 0 stress 0.1950753
## Run 1 stress 0.212001
## Run 2 stress 0.200252
## Run 3 stress 0.1996066
## Run 4 stress 0.1976846
## Run 5 stress 0.1995561
## Run 6 stress 0.2033033
## Run 7 stress 0.1950901
## ... Procrustes: rmse 0.01680846 max resid 0.1909933
## Run 8 stress 0.196432
## Run 9 stress 0.1968841
## Run 10 stress 0.199347
## Run 11 stress 0.1949618
## ... New best solution
## ... Procrustes: rmse 0.006644276 max resid 0.1063006
## Run 12 stress 0.2136834
## Run 13 stress 0.1962465
## Run 14 stress 0.1961008
## Run 15 stress 0.1969792
## Run 16 stress 0.1991054
## Run 17 stress 0.2014105
## Run 18 stress 0.1950813
## ... Procrustes: rmse 0.01574845 max resid 0.1912546
## Run 19 stress 0.1949615
## ... New best solution
## ... Procrustes: rmse 0.0005197987 max resid 0.006363654
## ... Similar to previous best
## Run 20 stress 0.1961116
## *** Solution reached
## [1] "unifrac_raw_PCoA"
## [1] "unifrac_raw_UMAP"
## [1] "wunifrac_raw_NMDS"
## Run 0 stress 0.161153
## Run 1 stress 0.1818377
## Run 2 stress 0.1772397
## Run 3 stress 0.1875217
## Run 4 stress 0.1737821
## Run 5 stress 0.195064
## Run 6 stress 0.1781208
## Run 7 stress 0.1799434
## Run 8 stress 0.1667898
## Run 9 stress 0.176366
## Run 10 stress 0.1945547
## Run 11 stress 0.1856628
## Run 12 stress 0.1836673
## Run 13 stress 0.168293
## Run 14 stress 0.1641538
## Run 15 stress 0.1998411
## Run 16 stress 0.1982581
## Run 17 stress 0.1795297
## Run 18 stress 0.1878017
## Run 19 stress 0.1669929
## Run 20 stress 0.1635855
## *** No convergence -- monoMDS stopping criteria:
## 13: stress ratio > sratmax
## 7: scale factor of the gradient < sfgrmin
## [1] "wunifrac_raw_PCoA"
## [1] "wunifrac_raw_UMAP"
## [1] "euclidean_vst_NMDS"
## Run 0 stress 0.1609013
## Run 1 stress 0.1626088
## Run 2 stress 0.1648398
## Run 3 stress 0.1649715
## Run 4 stress 0.164929
## Run 5 stress 0.1619001
## Run 6 stress 0.1635389
## Run 7 stress 0.1676298
## Run 8 stress 0.1636399
## Run 9 stress 0.1642341
## Run 10 stress 0.165878
## Run 11 stress 0.163316
## Run 12 stress 0.1691651
## Run 13 stress 0.1622519
## Run 14 stress 0.1628481
## Run 15 stress 0.1651573
## Run 16 stress 0.1664558
## Run 17 stress 0.1643859
## Run 18 stress 0.1637193
## Run 19 stress 0.1628497
## Run 20 stress 0.1603392
## ... New best solution
## ... Procrustes: rmse 0.02521003 max resid 0.1854854
## *** No convergence -- monoMDS stopping criteria:
## 18: no. of iterations >= maxit
## 2: stress ratio > sratmax
## [1] "euclidean_vst_PCoA"
## [1] "euclidean_vst_UMAP"
## [1] "bray_vst_NMDS"
## Run 0 stress 0.2277878
## Run 1 stress 0.2278856
## ... Procrustes: rmse 0.007424575 max resid 0.0859814
## Run 2 stress 0.2278728
## ... Procrustes: rmse 0.009270762 max resid 0.0848336
## Run 3 stress 0.2279974
## ... Procrustes: rmse 0.01041228 max resid 0.1189395
## Run 4 stress 0.2280305
## ... Procrustes: rmse 0.00752595 max resid 0.08207802
## Run 5 stress 0.2289447
## Run 6 stress 0.2324784
## Run 7 stress 0.2278391
## ... Procrustes: rmse 0.0116298 max resid 0.08629948
## Run 8 stress 0.2273921
## ... New best solution
## ... Procrustes: rmse 0.01005542 max resid 0.08419402
## Run 9 stress 0.2290124
## Run 10 stress 0.2318222
## Run 11 stress 0.2305865
## Run 12 stress 0.2275014
## ... Procrustes: rmse 0.009502425 max resid 0.0848474
## Run 13 stress 0.2285419
## Run 14 stress 0.2305028
## Run 15 stress 0.417808
## Run 16 stress 0.2295759
## Run 17 stress 0.2280828
## Run 18 stress 0.2303905
## Run 19 stress 0.2288548
## Run 20 stress 0.227982
## *** No convergence -- monoMDS stopping criteria:
## 1: no. of iterations >= maxit
## 19: stress ratio > sratmax
## [1] "bray_vst_PCoA"
## [1] "bray_vst_UMAP"
## [1] "jaccard_vst_NMDS"
## Run 0 stress 0.227261
## Run 1 stress 0.2328641
## Run 2 stress 0.2273415
## ... Procrustes: rmse 0.005575891 max resid 0.08523176
## Run 3 stress 0.2290474
## Run 4 stress 0.2273478
## ... Procrustes: rmse 0.005451431 max resid 0.08526879
## Run 5 stress 0.2305623
## Run 6 stress 0.2317595
## Run 7 stress 0.227494
## ... Procrustes: rmse 0.003573452 max resid 0.0507231
## Run 8 stress 0.2304353
## Run 9 stress 0.2280415
## Run 10 stress 0.2281836
## Run 11 stress 0.2309814
## Run 12 stress 0.2293575
## Run 13 stress 0.2272952
## ... Procrustes: rmse 0.006717131 max resid 0.083551
## Run 14 stress 0.2324858
## Run 15 stress 0.231942
## Run 16 stress 0.2292307
## Run 17 stress 0.2309612
## Run 18 stress 0.2306358
## Run 19 stress 0.2300717
## Run 20 stress 0.227541
## ... Procrustes: rmse 0.01020147 max resid 0.08705875
## *** No convergence -- monoMDS stopping criteria:
## 1: no. of iterations >= maxit
## 19: stress ratio > sratmax
## [1] "jaccard_vst_PCoA"
## [1] "jaccard_vst_UMAP"
## [1] "unifrac_vst_NMDS"
## Run 0 stress 1.886457e-16
## Run 1 stress 0
## ... New best solution
## ... Procrustes: rmse 0.05733635 max resid 0.869531
## Run 2 stress 0
## ... Procrustes: rmse 0.06038091 max resid 0.2822453
## Run 3 stress 0
## ... Procrustes: rmse 0.06056925 max resid 0.2859596
## Run 4 stress 0
## ... Procrustes: rmse 0.06056342 max resid 0.2914502
## Run 5 stress 0
## ... Procrustes: rmse 0.06068206 max resid 0.3025266
## Run 6 stress 0
## ... Procrustes: rmse 0.06090321 max resid 0.3274585
## Run 7 stress 0
## ... Procrustes: rmse 0.06035597 max resid 0.2793788
## Run 8 stress 0
## ... Procrustes: rmse 0.06047073 max resid 0.300135
## Run 9 stress 0
## ... Procrustes: rmse 0.06069685 max resid 0.2949206
## Run 10 stress 0
## ... Procrustes: rmse 0.06040829 max resid 0.2830932
## Run 11 stress 0
## ... Procrustes: rmse 0.0601835 max resid 0.2831408
## Run 12 stress 0
## ... Procrustes: rmse 0.06082349 max resid 0.3174186
## Run 13 stress 0
## ... Procrustes: rmse 0.06002378 max resid 0.2948504
## Run 14 stress 0
## ... Procrustes: rmse 0.05956902 max resid 0.2573516
## Run 15 stress 0
## ... Procrustes: rmse 0.05996597 max resid 0.2811964
## Run 16 stress 0
## ... Procrustes: rmse 0.06054728 max resid 0.2931258
## Run 17 stress 0
## ... Procrustes: rmse 0.06052762 max resid 0.3052316
## Run 18 stress 0
## ... Procrustes: rmse 0.06018397 max resid 0.27517
## Run 19 stress 0
## ... Procrustes: rmse 0.06056692 max resid 0.2856627
## Run 20 stress 0
## ... Procrustes: rmse 0.06072305 max resid 0.3204826
## *** No convergence -- monoMDS stopping criteria:
## 20: stress < smin
## [1] "unifrac_vst_PCoA"
## [1] "unifrac_vst_UMAP"
## [1] "wunifrac_vst_NMDS"
## Run 0 stress 0.1823021
## Run 1 stress 0.417799
## Run 2 stress 0.1823015
## ... New best solution
## ... Procrustes: rmse 0.0005246753 max resid 0.007854011
## ... Similar to previous best
## Run 3 stress 0.1858372
## Run 4 stress 0.1823021
## ... Procrustes: rmse 0.0004942914 max resid 0.007216553
## ... Similar to previous best
## Run 5 stress 0.1858381
## Run 6 stress 0.1866836
## Run 7 stress 0.1861078
## Run 8 stress 0.1823008
## ... New best solution
## ... Procrustes: rmse 0.0002517982 max resid 0.003748713
## ... Similar to previous best
## Run 9 stress 0.1865072
## Run 10 stress 0.1823024
## ... Procrustes: rmse 0.000312671 max resid 0.003092625
## ... Similar to previous best
## Run 11 stress 0.1866844
## Run 12 stress 0.4176636
## Run 13 stress 0.1823022
## ... Procrustes: rmse 0.0003406208 max resid 0.003818625
## ... Similar to previous best
## Run 14 stress 0.1858385
## Run 15 stress 0.1866838
## Run 16 stress 0.1823019
## ... Procrustes: rmse 0.0002926221 max resid 0.003797772
## ... Similar to previous best
## Run 17 stress 0.1823026
## ... Procrustes: rmse 0.0003126271 max resid 0.004146743
## ... Similar to previous best
## Run 18 stress 0.1858369
## Run 19 stress 0.1860942
## Run 20 stress 0.1827007
## ... Procrustes: rmse 0.004226909 max resid 0.06786374
## *** Solution reached
## [1] "wunifrac_vst_PCoA"
## [1] "wunifrac_vst_UMAP"
## [1] "euclidean_frac_NMDS"
## Run 0 stress 0.1204247
## Run 1 stress 0.1206017
## ... Procrustes: rmse 0.0400446 max resid 0.3173344
## Run 2 stress 0.1230261
## Run 3 stress 0.1254943
## Run 4 stress 0.1185844
## ... New best solution
## ... Procrustes: rmse 0.04744261 max resid 0.3797871
## Run 5 stress 0.121037
## Run 6 stress 0.125874
## Run 7 stress 0.1229305
## Run 8 stress 0.1258877
## Run 9 stress 0.1206976
## Run 10 stress 0.1187527
## ... Procrustes: rmse 0.04829711 max resid 0.4092316
## Run 11 stress 0.1241809
## Run 12 stress 0.1245341
## Run 13 stress 0.1308419
## Run 14 stress 0.1216227
## Run 15 stress 0.1219218
## Run 16 stress 0.1216191
## Run 17 stress 0.1250719
## Run 18 stress 0.1204038
## Run 19 stress 0.1221799
## Run 20 stress 0.1308524
## *** No convergence -- monoMDS stopping criteria:
## 19: no. of iterations >= maxit
## 1: stress ratio > sratmax
## [1] "euclidean_frac_PCoA"
## [1] "euclidean_frac_UMAP"
## [1] "bray_frac_NMDS"
## Run 0 stress 0.2818453
## Run 1 stress 0.2839507
## Run 2 stress 0.2833467
## Run 3 stress 0.2814406
## ... New best solution
## ... Procrustes: rmse 0.01850348 max resid 0.1865741
## Run 4 stress 0.2843496
## Run 5 stress 0.2836991
## Run 6 stress 0.2816992
## ... Procrustes: rmse 0.0195696 max resid 0.1104788
## Run 7 stress 0.2827069
## Run 8 stress 0.2930529
## Run 9 stress 0.283718
## Run 10 stress 0.2852771
## Run 11 stress 0.2856869
## Run 12 stress 0.2837538
## Run 13 stress 0.2819936
## Run 14 stress 0.2857617
## Run 15 stress 0.2827811
## Run 16 stress 0.2871906
## Run 17 stress 0.2831708
## Run 18 stress 0.2845995
## Run 19 stress 0.2855095
## Run 20 stress 0.2811038
## ... New best solution
## ... Procrustes: rmse 0.01760839 max resid 0.112224
## *** No convergence -- monoMDS stopping criteria:
## 4: no. of iterations >= maxit
## 16: stress ratio > sratmax
## [1] "bray_frac_PCoA"
## [1] "bray_frac_UMAP"
## [1] "jaccard_frac_NMDS"
## Run 0 stress 0.2825359
## Run 1 stress 0.295579
## Run 2 stress 0.2864782
## Run 3 stress 0.2829263
## ... Procrustes: rmse 0.02673403 max resid 0.1903155
## Run 4 stress 0.284254
## Run 5 stress 0.2838242
## Run 6 stress 0.2858169
## Run 7 stress 0.2886027
## Run 8 stress 0.2842911
## Run 9 stress 0.2869574
## Run 10 stress 0.2855451
## Run 11 stress 0.2858612
## Run 12 stress 0.2821634
## ... New best solution
## ... Procrustes: rmse 0.02578425 max resid 0.1844668
## Run 13 stress 0.2837911
## Run 14 stress 0.2868005
## Run 15 stress 0.2823502
## ... Procrustes: rmse 0.02333661 max resid 0.1862586
## Run 16 stress 0.282314
## ... Procrustes: rmse 0.02790461 max resid 0.130169
## Run 17 stress 0.2823868
## ... Procrustes: rmse 0.02675134 max resid 0.1843313
## Run 18 stress 0.2838576
## Run 19 stress 0.28254
## ... Procrustes: rmse 0.0246459 max resid 0.1156521
## Run 20 stress 0.2843962
## *** No convergence -- monoMDS stopping criteria:
## 7: no. of iterations >= maxit
## 13: stress ratio > sratmax
## [1] "jaccard_frac_PCoA"
## [1] "jaccard_frac_UMAP"
## [1] "unifrac_frac_NMDS"
## Run 0 stress 0.1950753
## Run 1 stress 0.1963941
## Run 2 stress 0.1960616
## Run 3 stress 0.1969395
## Run 4 stress 0.1982781
## Run 5 stress 0.1955889
## Run 6 stress 0.1969686
## Run 7 stress 0.2170798
## Run 8 stress 0.1947381
## ... New best solution
## ... Procrustes: rmse 0.009979722 max resid 0.09647186
## Run 9 stress 0.1977129
## Run 10 stress 0.1988329
## Run 11 stress 0.1984302
## Run 12 stress 0.1996111
## Run 13 stress 0.19627
## Run 14 stress 0.1968442
## Run 15 stress 0.1967855
## Run 16 stress 0.197732
## Run 17 stress 0.1999234
## Run 18 stress 0.2007899
## Run 19 stress 0.1968475
## Run 20 stress 0.1978205
## *** No convergence -- monoMDS stopping criteria:
## 1: no. of iterations >= maxit
## 16: stress ratio > sratmax
## 3: scale factor of the gradient < sfgrmin
## [1] "unifrac_frac_PCoA"
## [1] "unifrac_frac_UMAP"
## [1] "wunifrac_frac_NMDS"
## Run 0 stress 0.1693536
## Run 1 stress 0.1822301
## Run 2 stress 0.1862589
## Run 3 stress 0.1724381
## Run 4 stress 0.1700249
## Run 5 stress 0.1971008
## Run 6 stress 0.1890689
## Run 7 stress 0.1929433
## Run 8 stress 0.2079395
## Run 9 stress 0.1728318
## Run 10 stress 0.2044746
## Run 11 stress 0.1960163
## Run 12 stress 0.1776627
## Run 13 stress 0.1936799
## Run 14 stress 0.1741476
## Run 15 stress 0.1897843
## Run 16 stress 0.1836593
## Run 17 stress 0.1833076
## Run 18 stress 0.1710554
## Run 19 stress 0.1721027
## Run 20 stress 0.2050091
## *** No convergence -- monoMDS stopping criteria:
## 12: stress ratio > sratmax
## 8: scale factor of the gradient < sfgrmin
## [1] "wunifrac_frac_PCoA"
## [1] "wunifrac_frac_UMAP"
From the plots it looks like the vst normalization in the UMAP clustering for euclidean distances gives best results for separation of Hydra from different locations directly. If one settles on one of the mostly used algorithms it looks like vst-bray and vst-jaccard give reliable results (at least UMAP is able to cluster them quite well). The data has many dimensions, which makes it difficult for NMDS an PCoA to reveal common patterns in the samples. Surprisingly the unifrac and wunifrac distances look rather scattered and give relatively bad clustering, which might imply that the inferred fast-tree ML-tree is no really representing the phylogenetic distance of the bacterial species, or that that the phylogenetic distance is not a good descriptor for the samples.
Lets see, whether the ordination is also representative for Species and Reproductive Mode differences.
sel <- paste(c("euclidean","bray","jaccard"),"vst", sep ="_")
plots2 <- list()
for(nm in sel){
ordi <- umap(as.matrix(mtrxs[[nm]]), input="dist")
plots2[[paste(nm,"RepMod",sep="_")]] <- ggplot(as.data.frame(ordi$layout), aes(x=V1,y=V2, color = meta$ReproductiveModeSimple))+
ggtitle(nm) +
labs(color = "RepMode")+
geom_point()
plots2[[paste(nm,"Species",sep="_")]] <- ggplot(as.data.frame(ordi$layout), aes(x=V1,y=V2, color = meta$Species))+
ggtitle(nm) +
labs(color = "Species")+
geom_point()
plots2[[paste(nm,"PopID",sep="_")]] <- ggplot(as.data.frame(ordi$layout), aes(x=V1,y=V2, color = meta$PopID))+
ggtitle(nm) +
labs(color = "PopID")+
geom_point()
}
pp_UMAP <- cowplot::plot_grid(plotlist=plots2, ncol = 3)
pp_UMAP
Knowing these things, I will calculate an adonis on the euclidean, bray and jaccard distances for the vst normalized data and check for significant differences in sampling location (= Population ID [PopID]), reproductive mode (ReproductiveMode) and Species (Species).
sel <- paste(c("euclidean","bray","jaccard"),"vst", sep ="_")
stats <- list()
for(n in sel){
stats[[n]] <- adonis(formula = mtrxs[[n]] ~meta$Species+meta$ReproductiveModeSimple+ meta$PopID, perm = 200)
}
print("Adonis testing for PopID")## [1] "Adonis testing for PopID"print(stats[sel])## $euclidean_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$Species + meta$ReproductiveModeSimple + meta$PopID, permutations = 200)
##
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$Species 2 19538 9769.2 4.4897 0.02572 0.004975 **
## meta$ReproductiveModeSimple 2 8230 4115.1 1.8912 0.01083 0.004975 **
## meta$PopID 20 209706 10485.3 4.8188 0.27604 0.004975 **
## Residuals 240 522223 2175.9 0.68741
## Total 264 759698 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $bray_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$Species + meta$ReproductiveModeSimple + meta$PopID, permutations = 200)
##
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$Species 2 3.029 1.51455 7.9707 0.03920 0.004975 **
## meta$ReproductiveModeSimple 2 0.918 0.45886 2.4149 0.01188 0.004975 **
## meta$PopID 20 27.732 1.38660 7.2974 0.35884 0.004975 **
## Residuals 240 45.603 0.19001 0.59009
## Total 264 77.282 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $jaccard_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$Species + meta$ReproductiveModeSimple + meta$PopID, permutations = 200)
##
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$Species 2 3.019 1.50965 5.3374 0.03081 0.004975 **
## meta$ReproductiveModeSimple 2 1.030 0.51505 1.8210 0.01051 0.004975 **
## meta$PopID 20 26.052 1.30260 4.6053 0.26588 0.004975 **
## Residuals 240 67.883 0.28285 0.69279
## Total 264 97.984 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1stats0.5 <- list()
for(n in sel){
stats0.5[[n]] <- adonis(formula = mtrxs[[n]] ~meta$PopID, perm = 200)
}
print("Adonis testing for PopID ONLY")## [1] "Adonis testing for PopID ONLY"print(stats0.5[sel])## $euclidean_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$PopID, permutations = 200)
##
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$PopID 20 223721 11186.1 5.0924 0.29449 0.004975 **
## Residuals 244 535977 2196.6 0.70551
## Total 264 759698 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $bray_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$PopID, permutations = 200)
##
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$PopID 20 29.950 1.49749 7.7196 0.38754 0.004975 **
## Residuals 244 47.332 0.19398 0.61246
## Total 264 77.282 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $jaccard_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$PopID, permutations = 200)
##
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$PopID 20 28.012 1.40060 4.884 0.28588 0.004975 **
## Residuals 244 69.972 0.28677 0.71412
## Total 264 97.984 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1stats2 <- list()
for(n in sel){
stats2[[n]] <- adonis(formula = mtrxs[[n]] ~ meta$Species, strata= meta$PopID, perm =200)
}
print("Adonis testing for Species and Reproductive Mode")## [1] "Adonis testing for Species and Reproductive Mode"print(stats2)## $euclidean_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$Species, permutations = 200, strata = meta$PopID)
##
## Blocks: strata
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$Species 2 19538 9769.2 3.4581 0.02572 0.004975 **
## Residuals 262 740160 2825.0 0.97428
## Total 264 759698 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $bray_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$Species, permutations = 200, strata = meta$PopID)
##
## Blocks: strata
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$Species 2 3.029 1.51455 5.3441 0.0392 0.004975 **
## Residuals 262 74.253 0.28341 0.9608
## Total 264 77.282 1.0000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $jaccard_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta$Species, permutations = 200, strata = meta$PopID)
##
## Blocks: strata
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta$Species 2 3.019 1.50965 4.165 0.03081 0.004975 **
## Residuals 262 94.965 0.36246 0.96919
## Total 264 97.984 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# select only oligactis for species test in adonis
sel2 <- rownames(meta[meta$Species == "OLI",])
stats3 <- list()
for(n in sel){
stats3[[n]] <- adonis(formula = as.matrix(mtrxs[[n]])[sel2,sel2] ~ meta[sel2,"ReproductiveModeSimple"], strata= meta[sel2,"PopID"], perm =200)
}
print("Adonis testing for Reproductive Mode")## [1] "Adonis testing for Reproductive Mode"print(stats3)## $euclidean_vst
##
## Call:
## adonis(formula = as.matrix(mtrxs[[n]])[sel2, sel2] ~ meta[sel2, "ReproductiveModeSimple"], permutations = 200, strata = meta[sel2, "PopID"])
##
## Blocks: strata
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2
## meta[sel2, "ReproductiveModeSimple"] 2 8536 4268.0 1.376 0.01318
## Residuals 206 638971 3101.8 0.98682
## Total 208 647507 1.00000
## Pr(>F)
## meta[sel2, "ReproductiveModeSimple"] 0.4826
## Residuals
## Total
##
## $bray_vst
##
## Call:
## adonis(formula = as.matrix(mtrxs[[n]])[sel2, sel2] ~ meta[sel2, "ReproductiveModeSimple"], permutations = 200, strata = meta[sel2, "PopID"])
##
## Blocks: strata
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2
## meta[sel2, "ReproductiveModeSimple"] 2 0.973 0.48635 1.6038 0.01533
## Residuals 206 62.468 0.30324 0.98467
## Total 208 63.441 1.00000
## Pr(>F)
## meta[sel2, "ReproductiveModeSimple"] 0.408
## Residuals
## Total
##
## $jaccard_vst
##
## Call:
## adonis(formula = as.matrix(mtrxs[[n]])[sel2, sel2] ~ meta[sel2, "ReproductiveModeSimple"], permutations = 200, strata = meta[sel2, "PopID"])
##
## Blocks: strata
## Permutation: free
## Number of permutations: 200
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2
## meta[sel2, "ReproductiveModeSimple"] 2 1.062 0.53086 1.4001 0.01341
## Residuals 206 78.109 0.37917 0.98659
## Total 208 79.171 1.00000
## Pr(>F)
## meta[sel2, "ReproductiveModeSimple"] 0.4478
## Residuals
## TotalTo finally answer the main question from the beginning: 1.) There is a consistent and rather strong signature in the data which is indicative for the sampling location, hence bacterial communities in Hydra differ if they live in different locations. 2.) There is a consistent difference in the bacterial communities from different Hydra species (this is what was expected). 3.) There is a less strong, but detectable change in the microbial community if the animals have different reproductive modes. This is only detectable if the model design is stratified by the sampling location of the 16S data. However, it indicates, that change in the reproductive mode is associated with the differences in the microbiota.
We obtained new meta data for the lakes where the animals were sampled and Sebastian suggested to test for differences in the microbiota in dependence of the water body. Does it make a difference if Hydra lives in a lake or a flowing (river) water body? Further we obtained the geocoordinates of the sampling sites. I will use the distance between geolocation and \(\beta\)-diversity to correlate these two frameworks.
metaLakes <- read.csv(snakemake@input[["MetaLakes"]])
# join the meta data
meta2 <- merge(meta, metaLakes, by.x = "PopID", by.y = "SiteID", all.x=TRUE)
stats3 <- list()
for(n in sel){
stats3[[n]] <- adonis(formula = mtrxs[[n]] ~ meta2$waterbody, perm =500)
}
stats3## $euclidean_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta2$waterbody, permutations = 500)
##
## Permutation: free
## Number of permutations: 500
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta2$waterbody 1 9838 9838.3 3.4506 0.01295 0.001996 **
## Residuals 263 749860 2851.2 0.98705
## Total 264 759698 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $bray_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta2$waterbody, permutations = 500)
##
## Permutation: free
## Number of permutations: 500
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta2$waterbody 1 1.110 1.10971 3.8315 0.01436 0.001996 **
## Residuals 263 76.172 0.28963 0.98564
## Total 264 77.282 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $jaccard_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ meta2$waterbody, permutations = 500)
##
## Permutation: free
## Number of permutations: 500
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## meta2$waterbody 1 1.082 1.08223 2.9373 0.01104 0.001996 **
## Residuals 263 96.902 0.36845 0.98896
## Total 264 97.984 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# The waterbody is a factor which changes the microbiota!
stats4 <- list()
for(n in sel){
stats4[[n]] <- adonis(formula = mtrxs[[n]] ~ factor(meta2$NutrientLoad), perm =500)
}
stats4## $euclidean_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ factor(meta2$NutrientLoad), permutations = 500)
##
## Permutation: free
## Number of permutations: 500
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## factor(meta2$NutrientLoad) 2 12479 6239.4 2.1877 0.01643 0.001996 **
## Residuals 262 747219 2852.0 0.98357
## Total 264 759698 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $bray_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ factor(meta2$NutrientLoad), permutations = 500)
##
## Permutation: free
## Number of permutations: 500
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## factor(meta2$NutrientLoad) 2 1.694 0.84709 2.9362 0.02192 0.001996 **
## Residuals 262 75.588 0.28850 0.97808
## Total 264 77.282 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $jaccard_vst
##
## Call:
## adonis(formula = mtrxs[[n]] ~ factor(meta2$NutrientLoad), permutations = 500)
##
## Permutation: free
## Number of permutations: 500
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## factor(meta2$NutrientLoad) 2 1.697 0.84834 2.3083 0.01732 0.001996 **
## Residuals 262 96.288 0.36751 0.98268
## Total 264 97.984 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# calculate the geographic distance between sampling points
geoDist <- distm(metaLakes[,c("lon","lat")])
colnames(geoDist) <- metaLakes[["SiteID"]]
rownames(geoDist) <- metaLakes[["SiteID"]]
# set the upper triangle to NA, so its easily filtered afterwards
geoDist[upper.tri(geoDist, diag=FALSE)] <- NA
# create a data frame which associates distance of microbiome to geographic distance
geoDistLong <- data.table(reshape2::melt(geoDist))
geoDistLong <- geoDistLong[!is.na(value),]
mtrxsLong <- list()
for(mat in sel){
bacDist <- as.matrix(mtrxs[[mat]])
bacDist[upper.tri(bacDist, diag=TRUE)] <- NA
bacDistLong <- data.table(reshape2::melt(bacDist))
bacDistLong <- bacDistLong[!is.na(value),]
# add the PopID information to associate it back to the geographic distance
bacDistLong <- cbind(bacDistLong,
Var1PopID = meta[bacDistLong[,Var1], "PopID"],
Var2PopID = meta[bacDistLong[,Var2], "PopID"],
GeoDist = -1)
mtrxsLong[[mat]] <- bacDistLong
}
# merge both tables
for(i in 1:nrow(geoDistLong)){
for(mat in names(mtrxsLong)){
sel2 <- ((mtrxsLong[[mat]][,Var1PopID] %in% geoDistLong[i,Var1] &
mtrxsLong[[mat]][,Var2PopID] %in% geoDistLong[i,Var2]) |
(mtrxsLong[[mat]][,Var2PopID] %in% geoDistLong[i,Var1] &
mtrxsLong[[mat]][,Var1PopID] %in% geoDistLong[i,Var2]))
mtrxsLong[[mat]][sel2, "GeoDist"] <- geoDistLong[i,value]
}
}
# create simple dotplots
geoDistPlots <- list()
geoDistModels <- list()
for(mat in names(mtrxsLong)){
nm <- strsplit(mat, split="_")[[1]]
p <- ggplot(mtrxsLong[[mat]], aes(y=value, x=GeoDist)) +
geom_point(alpha = 0.1) +
xlab("Geographic distance") +
ylab(paste(nm[1], "distance")) +
ggtitle(paste("Distance based on", nm[2]))
geoDistPlots[[mat]] <- p
print(cor.test(method ="pearson", mtrxsLong[[mat]][,value], mtrxsLong[[mat]][,GeoDist]))
geoDistModels[[mat]] <- lm(data=mtrxsLong[[mat]], value~GeoDist)
}##
## Pearson's product-moment correlation
##
## data: mtrxsLong[[mat]][, value] and mtrxsLong[[mat]][, GeoDist]
## t = -14.565, df = 34978, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.08804867 -0.06721599
## sample estimates:
## cor
## -0.07764081
##
##
## Pearson's product-moment correlation
##
## data: mtrxsLong[[mat]][, value] and mtrxsLong[[mat]][, GeoDist]
## t = 25.377, df = 34978, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.1241526 0.1447327
## sample estimates:
## cor
## 0.1344572
##
##
## Pearson's product-moment correlation
##
## data: mtrxsLong[[mat]][, value] and mtrxsLong[[mat]][, GeoDist]
## t = 27.216, df = 34978, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.1337289 0.1542533
## sample estimates:
## cor
## 0.1440066pp <- cowplot::plot_grid(plotlist = geoDistPlots, ncol = 1)
pp
# create a distance matrix for the distance in geo-location of the same dimensions like the beta-diversity to perform a mantel test
# 1. cast into wide format
geoDist_Samples <- dcast(mtrxsLong[[1]], Var1~Var2, Value.var = "GeoDist")
# 2. format as data.frame, correct rownames, which are missing in data.table
geoDist_Samples_names <- geoDist_Samples[, Var1]
geoDist_Samples <- as.data.frame(geoDist_Samples)
rownames(geoDist_Samples) <- geoDist_Samples_names
geoDist_Samples <- geoDist_Samples[,-1]
# 3. add an additional column and row, because I left out the diagonal by the first conversion from distance to long format
geoDist_Samples <- cbind(geoDist_Samples, added = NA)
colnames(geoDist_Samples)[ncol(geoDist_Samples)] <- geoDist_Samples_names[length(geoDist_Samples_names)]
geoDist_Samples_names <- colnames(geoDist_Samples)
geoDist_Samples <- rbind(added = NA, geoDist_Samples)
rownames(geoDist_Samples)[1] <- geoDist_Samples_names[1]
# 4. reformat into distance object
geoDist_Samples <- as.dist(geoDist_Samples)
# perform a mantel test on the distance matrices
for(mat in names(mtrxsLong)){
print(mat)
print(mantel(geoDist_Samples, mtrxs[[mat]]))
}## [1] "euclidean_vst"
##
## Mantel statistic based on Pearson's product-moment correlation
##
## Call:
## mantel(xdis = geoDist_Samples, ydis = mtrxs[[mat]])
##
## Mantel statistic r: 1
## Significance: 0.001
##
## Upper quantiles of permutations (null model):
## 90% 95% 97.5% 99%
## 0.0770 0.0981 0.1168 0.1451
## Permutation: free
## Number of permutations: 999
##
## [1] "bray_vst"
##
## Mantel statistic based on Pearson's product-moment correlation
##
## Call:
## mantel(xdis = geoDist_Samples, ydis = mtrxs[[mat]])
##
## Mantel statistic r: 0.7914
## Significance: 0.001
##
## Upper quantiles of permutations (null model):
## 90% 95% 97.5% 99%
## 0.0553 0.0727 0.0867 0.1004
## Permutation: free
## Number of permutations: 999
##
## [1] "jaccard_vst"
##
## Mantel statistic based on Pearson's product-moment correlation
##
## Call:
## mantel(xdis = geoDist_Samples, ydis = mtrxs[[mat]])
##
## Mantel statistic r: 0.7649
## Significance: 0.001
##
## Upper quantiles of permutations (null model):
## 90% 95% 97.5% 99%
## 0.0583 0.0726 0.0860 0.0959
## Permutation: free
## Number of permutations: 999So there is a high statistic significant association of geographic distance and \(\beta\)-diversity, however, the effect is marginal and probably only significant because of the shear amount of data in the set. Its a nice feature, but I would not make to much a deal out of it.
# The idea of Jacint was to associate the alpha diversity to the
alphaDiv <- estimate_richness(phy)
alphaNames <- colnames(alphaDiv)
alphaDiv <- merge(alphaDiv, meta, by=0, all.x=TRUE)
# transform the data to fit normality a little better
alphaDiv[["Simpson"]] <- alphaDiv[["Simpson"]]^2
alphaDiv[["Chao1"]] <- log(alphaDiv[["Chao1"]])
alphaPlots <- list()
for(nm in c("Simpson","Shannon","Chao1")){
for(group in c("PopID", "Species", "ReproductiveModeSimple")){
id <- paste(nm, group, sep ="_")
p <- ggplot(alphaDiv, aes_string(x = group, y= nm, fill = group))+
geom_boxplot() +
ylab(nm) +
xlab(group) +
theme(legend.position = "None")
alphaPlots[[id]] <- p
}
}
pp <- cowplot::plot_grid(plotlist = alphaPlots, rel_widths = rep(c(3,1,1),3))
print(pp)
There seems to be major differences in the alpha diversity of the different populations, but also in the different species. Lets fit some linear (mixed) models to the data and see if these differences are significant.
Lets start with the Simpson index:
# lets fit a linear mixed model to the data and check for population
simpFull <- lm(data=alphaDiv, Simpson ~ PopID+Species+ReproductiveModeSimple)
simpPopSpec <- lm(data=alphaDiv, Simpson ~ PopID+Species)
simpPopNull <- lm(data=alphaDiv, Simpson ~ PopID)
simpPopLMM <- lme(data=alphaDiv, Simpson ~ PopID, random = ~1|Species/ReproductiveModeSimple, method = "ML")
anova(simpPopLMM, simpFull)## Model df AIC BIC logLik Test L.Ratio p-value
## simpPopLMM 1 24 -22.26486 63.64866 35.13243
## simpFull 2 26 -32.50277 60.57021 42.25138 1 vs 2 14.23791 8e-04anova(simpPopLMM, simpPopNull)## Model df AIC BIC logLik Test L.Ratio p-value
## simpPopLMM 1 24 -22.26486 63.64866 35.13243
## simpPopNull 2 22 -15.98438 62.76967 29.99219 1 vs 2 10.28047 0.0059# the full model is still better, test the difference to the pop+species model
anova(simpPopNull,simpPopSpec,simpFull)## Analysis of Variance Table
##
## Model 1: Simpson ~ PopID
## Model 2: Simpson ~ PopID + Species
## Model 3: Simpson ~ PopID + Species + ReproductiveModeSimple
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 244 12.373
## 2 242 11.553 2 0.81973 8.7210 0.0002207 ***
## 3 240 11.279 2 0.27367 2.9115 0.0563182 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(simpPopNull,simpPopSpec,simpFull)## df AIC
## simpPopNull 22 -15.98438
## simpPopSpec 24 -30.14997
## simpFull 26 -32.50277BIC(simpPopNull,simpPopSpec,simpFull)## df BIC
## simpPopNull 22 62.76967
## simpPopSpec 24 55.76355
## simpFull 26 60.57021# the linear model with only Population and Species is best fitting the data
# lets check for species differences
simpSpecNull <- lm(data=alphaDiv, Simpson ~ Species)
simpSpecLMM1 <- lme(data=alphaDiv, Simpson ~ Species, random =~1|PopID, method = "ML")
simpSpecLMM2 <- lme(data=alphaDiv, Simpson ~ Species+ReproductiveModeSimple, random =~1|PopID, method = "ML")
anova(simpSpecLMM2, simpSpecLMM1)## Model df AIC BIC logLik Test L.Ratio p-value
## simpSpecLMM2 1 7 -5.466983 19.59113 9.733491
## simpSpecLMM1 2 5 -4.822704 13.07595 7.411352 1 vs 2 4.644279 0.0981anova(simpSpecLMM2, simpFull)## Model df AIC BIC logLik Test L.Ratio p-value
## simpSpecLMM2 1 7 -5.46698 19.59113 9.73349
## simpFull 2 26 -32.50277 60.57021 42.25138 1 vs 2 65.03578 <.0001anova(simpSpecLMM2, simpSpecNull)## Model df AIC BIC logLik Test L.Ratio p-value
## simpSpecLMM2 1 7 -5.46698 19.59113 9.733491
## simpSpecNull 2 4 36.93405 51.25297 -14.467024 1 vs 2 48.40103 <.0001anova(simpSpecNull, simpPopSpec, simpFull)## Analysis of Variance Table
##
## Model 1: Simpson ~ Species
## Model 2: Simpson ~ PopID + Species
## Model 3: Simpson ~ PopID + Species + ReproductiveModeSimple
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 262 17.306
## 2 242 11.553 20 5.7527 6.1203 4.551e-13 ***
## 3 240 11.279 2 0.2737 2.9115 0.05632 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(simpSpecNull,simpPopSpec,simpFull)## df AIC
## simpSpecNull 4 36.93405
## simpPopSpec 24 -30.14997
## simpFull 26 -32.50277BIC(simpSpecNull,simpPopSpec,simpFull)## df BIC
## simpSpecNull 4 51.25297
## simpPopSpec 24 55.76355
## simpFull 26 60.57021# again the population+Species simple linear model is best fitting the data
simpMod <- simpPopSpec
# ok lets repeat this procedure with the other diversity measures
shanFull <- lm(data=alphaDiv, Shannon ~ PopID+Species+ReproductiveModeSimple)
shanPopSpec <- lm(data=alphaDiv, Shannon ~ PopID+Species)
shanPopNull <- lm(data= alphaDiv, Shannon ~ PopID)
shanSpecNull <- lm(data= alphaDiv, Shannon ~ Species)
shanPopLMM <- lme(data=alphaDiv, Shannon ~ PopID, random = ~1|Species/ReproductiveModeSimple, method = "ML")
shanSpecLMM1 <- lme(data=alphaDiv, Shannon~ Species, random =~1|PopID, method = "ML")
shanSpecLMM2 <- lme(data=alphaDiv, Shannon~ Species+ReproductiveModeSimple, random =~1|PopID, method = "ML")
# find the best population model
anova(shanPopLMM, shanPopNull)## Model df AIC BIC logLik Test L.Ratio p-value
## shanPopLMM 1 24 634.9987 720.9122 -293.4993
## shanPopNull 2 22 638.2704 717.0244 -297.1352 1 vs 2 7.271693 0.0264anova(shanPopLMM, shanFull)## Model df AIC BIC logLik Test L.Ratio p-value
## shanPopLMM 1 24 634.9987 720.9122 -293.4993
## shanFull 2 26 625.7518 718.8248 -286.8759 1 vs 2 13.24687 0.0013anova(shanPopNull, shanPopSpec, shanFull)## Analysis of Variance Table
##
## Model 1: Shannon ~ PopID
## Model 2: Shannon ~ PopID + Species
## Model 3: Shannon ~ PopID + Species + ReproductiveModeSimple
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 244 146.11
## 2 242 139.26 2 6.8489 6.0777 0.002662 **
## 3 240 135.23 2 4.0375 3.5829 0.029292 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(shanPopNull, shanPopSpec, shanFull)## df AIC
## shanPopNull 22 638.2704
## shanPopSpec 24 629.5482
## shanFull 26 625.7518BIC(shanPopNull, shanPopSpec, shanFull)## df BIC
## shanPopNull 22 717.0244
## shanPopSpec 24 715.4617
## shanFull 26 718.8248# this one is a bit tricky, but I would go to the full model here, as log-likelihood and AIC are in concordance here
# species effect
anova(shanSpecLMM1, shanSpecLMM2)## Model df AIC BIC logLik Test L.Ratio p-value
## shanSpecLMM1 1 5 655.3448 673.2434 -322.6724
## shanSpecLMM2 2 7 653.1819 678.2400 -319.5910 1 vs 2 6.162849 0.0459anova(shanSpecLMM2, shanSpecNull)## Model df AIC BIC logLik Test L.Ratio p-value
## shanSpecLMM2 1 7 653.1819 678.2400 -319.5910
## shanSpecNull 2 4 707.2619 721.5808 -349.6309 1 vs 2 60.07994 <.0001anova(shanSpecLMM2, shanFull)## Model df AIC BIC logLik Test L.Ratio p-value
## shanSpecLMM2 1 7 653.1819 678.2400 -319.5910
## shanFull 2 26 625.7518 718.8248 -286.8759 1 vs 2 65.43015 <.0001anova(shanSpecNull, shanPopSpec, shanFull)## Analysis of Variance Table
##
## Model 1: Shannon ~ Species
## Model 2: Shannon ~ PopID + Species
## Model 3: Shannon ~ PopID + Species + ReproductiveModeSimple
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 262 217.15
## 2 242 139.26 20 77.883 6.9113 5.385e-15 ***
## 3 240 135.23 2 4.038 3.5829 0.02929 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(shanSpecNull, shanPopSpec, shanFull)## df AIC
## shanSpecNull 4 707.2619
## shanPopSpec 24 629.5482
## shanFull 26 625.7518BIC(shanSpecNull, shanPopSpec, shanFull)## df BIC
## shanSpecNull 4 721.5808
## shanPopSpec 24 715.4617
## shanFull 26 718.8248# again the full model seems to be the best, at least concord in log-likelihood and AIC
shanMod <- shanFull
chaoFull <- lm(data=alphaDiv, Chao1 ~ PopID+Species+ReproductiveModeSimple)
chaoPopNull <- lm(data=alphaDiv, Chao1 ~ PopID)
chaoPopSpec <- lm(data=alphaDiv, Chao1 ~ PopID+Species)
chaoSpecNull<- lm(data=alphaDiv, Chao1 ~ Species)
chaoPopLMM <- lme(data=alphaDiv, Chao1 ~ PopID, random = ~1|Species/ReproductiveModeSimple, method = "ML")
chaoSpecLMM1 <- lme(data=alphaDiv, Chao1 ~ Species, random =~1|PopID, method = "ML")
chaoSpecLMM2 <- lme(data=alphaDiv, Chao1 ~ Species+ReproductiveModeSimple, random =~1|PopID, method = "ML")
# population
anova(chaoPopLMM, chaoPopNull)## Model df AIC BIC logLik Test L.Ratio p-value
## chaoPopLMM 1 24 447.1783 533.0919 -199.5892
## chaoPopNull 2 22 451.5532 530.3073 -203.7766 1 vs 2 8.374895 0.0152anova(chaoPopLMM, chaoFull)## Model df AIC BIC logLik Test L.Ratio p-value
## chaoPopLMM 1 24 447.1783 533.0919 -199.5892
## chaoFull 2 26 438.2494 531.3224 -193.1247 1 vs 2 12.92893 0.0016anova(chaoPopNull, chaoPopSpec, chaoFull)## Analysis of Variance Table
##
## Model 1: Chao1 ~ PopID
## Model 2: Chao1 ~ PopID + Species
## Model 3: Chao1 ~ PopID + Species + ReproductiveModeSimple
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 244 72.225
## 2 242 67.124 2 5.1012 9.1851 0.0001433 ***
## 3 240 66.646 2 0.4778 0.8603 0.4243201
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(chaoPopNull, chaoPopSpec, chaoFull)## df AIC
## chaoPopNull 22 451.5532
## chaoPopSpec 24 436.1425
## chaoFull 26 438.2494BIC(chaoPopNull, chaoPopSpec, chaoFull)## df BIC
## chaoPopNull 22 530.3073
## chaoPopSpec 24 522.0560
## chaoFull 26 531.3224# the linear model with population and species is fitting the data best
chaoMod <- chaoPopSpec
# species
anova(chaoSpecLMM1,chaoSpecLMM2 , chaoSpecNull)## Model df AIC BIC logLik Test L.Ratio p-value
## chaoSpecLMM1 1 5 477.4263 495.3250 -233.7132
## chaoSpecLMM2 2 7 479.7240 504.7821 -232.8620 1 vs 2 1.70232 0.4269
## chaoSpecNull 3 4 580.3695 594.6884 -286.1848 2 vs 3 106.64548 <.0001anova(chaoSpecLMM1,chaoSpecLMM2 , chaoFull)## Model df AIC BIC logLik Test L.Ratio p-value
## chaoSpecLMM1 1 5 477.4263 495.3250 -233.7132
## chaoSpecLMM2 2 7 479.7240 504.7821 -232.8620 1 vs 2 1.70232 0.4269
## chaoFull 3 26 438.2494 531.3224 -193.1247 2 vs 3 79.47462 <.0001anova(chaoSpecNull, chaoPopSpec, chaoFull)## Analysis of Variance Table
##
## Model 1: Chao1 ~ Species
## Model 2: Chao1 ~ PopID + Species
## Model 3: Chao1 ~ PopID + Species + ReproductiveModeSimple
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 262 134.523
## 2 242 67.124 20 67.399 12.1356 <2e-16 ***
## 3 240 66.646 2 0.478 0.8603 0.4243
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(chaoSpecNull, chaoPopSpec, chaoFull)## df AIC
## chaoSpecNull 4 580.3695
## chaoPopSpec 24 436.1425
## chaoFull 26 438.2494BIC(chaoSpecNull, chaoPopSpec, chaoFull)## df BIC
## chaoSpecNull 4 594.6884
## chaoPopSpec 24 522.0560
## chaoFull 26 531.3224OK after fitting the data to a (rather simple) model lets see if the assumptions hold an plot some diagnostics.
simp_pp <- plot_grid(~plot(simpMod, which =1),
~plot(simpMod, which =2),
~plot(simpMod,which=4),
~plot(simpMod,which=5))
shan_pp <- plot_grid(~plot(shanMod, which =1),
~plot(shanMod, which =2),
~plot(shanMod,which=4),
~plot(shanMod,which=5))
chao_pp <- plot_grid(~plot(chaoMod, which =1),
~plot(chaoMod, which =2),
~plot(chaoMod,which=4),
~plot(chaoMod,which=5))
simp_pp
shan_pp
chao_pp
The models look quite okay, lets look on the coefficients.
summary(simpMod)##
## Call:
## lm(formula = Simpson ~ PopID + Species, data = alphaDiv)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.63243 -0.12389 0.01434 0.13754 0.51192
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.18747 0.09809 1.911 0.057170 .
## PopIDM108 0.02174 0.09559 0.227 0.820279
## PopIDM109 0.27082 0.09783 2.768 0.006072 **
## PopIDM26 0.18345 0.10050 1.825 0.069189 .
## PopIDM28 0.20177 0.08539 2.363 0.018920 *
## PopIDM31 0.17151 0.10769 1.593 0.112574
## PopIDM34 0.19550 0.11293 1.731 0.084699 .
## PopIDM44 0.19147 0.08821 2.171 0.030928 *
## PopIDM67 -0.36262 0.09771 -3.711 0.000256 ***
## PopIDM70 -0.08486 0.10815 -0.785 0.433406
## PopIDM71 -0.09660 0.09203 -1.050 0.294922
## PopIDM72 0.03747 0.08408 0.446 0.656265
## PopIDM78 -0.18167 0.16931 -1.073 0.284339
## PopIDM79 0.07531 0.08931 0.843 0.399930
## PopIDM83 -0.03842 0.09090 -0.423 0.672861
## PopIDM84 0.20800 0.09203 2.260 0.024696 *
## PopIDM85 0.26310 0.08933 2.945 0.003541 **
## PopIDM86 0.55867 0.15980 3.496 0.000562 ***
## PopIDM89 0.13784 0.08809 1.565 0.118969
## PopIDM90 0.16183 0.08287 1.953 0.051991 .
## PopIDR12 -0.07562 0.09783 -0.773 0.440283
## SpeciesOLI 0.26641 0.07092 3.756 0.000216 ***
## SpeciesVUL 0.29167 0.07246 4.025 7.62e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2185 on 242 degrees of freedom
## Multiple R-squared: 0.3467, Adjusted R-squared: 0.2873
## F-statistic: 5.836 on 22 and 242 DF, p-value: 3.653e-13Anova(simpMod)## Anova Table (Type II tests)
##
## Response: Simpson
## Sum Sq Df F value Pr(>F)
## PopID 5.7527 20 6.0251 7.434e-13 ***
## Species 0.8197 2 8.5854 0.0002499 ***
## Residuals 11.5530 242
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1summary(shanMod)##
## Call:
## lm(formula = Shannon ~ PopID + Species + ReproductiveModeSimple,
## data = alphaDiv)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.46721 -0.47374 0.00618 0.40648 1.92327
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.4848 0.3661 1.324 0.186773
## PopIDM108 0.1324 0.3345 0.396 0.692450
## PopIDM109 1.1456 0.3361 3.408 0.000766 ***
## PopIDM26 1.1531 0.3530 3.266 0.001250 **
## PopIDM28 0.7970 0.3009 2.649 0.008607 **
## PopIDM31 0.6990 0.3738 1.870 0.062734 .
## PopIDM34 0.5082 0.3883 1.309 0.191765
## PopIDM44 0.7452 0.3109 2.397 0.017309 *
## PopIDM67 -1.1852 0.3445 -3.441 0.000684 ***
## PopIDM70 -0.1457 0.3733 -0.390 0.696635
## PopIDM71 -0.3625 0.3246 -1.117 0.265133
## PopIDM72 0.4837 0.2958 1.635 0.103288
## PopIDM78 -0.6623 0.5864 -1.129 0.259867
## PopIDM79 0.2035 0.3114 0.654 0.514057
## PopIDM83 0.1703 0.3223 0.528 0.597673
## PopIDM84 0.6658 0.3256 2.045 0.041942 *
## PopIDM85 1.4690 0.3124 4.703 4.33e-06 ***
## PopIDM86 1.4996 0.5514 2.720 0.007011 **
## PopIDM89 0.6154 0.3090 1.991 0.047575 *
## PopIDM90 0.4680 0.2883 1.623 0.105886
## PopIDR12 -0.2334 0.3363 -0.694 0.488338
## SpeciesOLI 1.1089 0.2634 4.210 3.62e-05 ***
## SpeciesVUL 1.0162 0.2750 3.695 0.000272 ***
## ReproductiveModeSimpleNR 0.1121 0.1493 0.751 0.453500
## ReproductiveModeSimpleSEX 0.3651 0.1364 2.677 0.007944 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7506 on 240 degrees of freedom
## Multiple R-squared: 0.3936, Adjusted R-squared: 0.3329
## F-statistic: 6.49 on 24 and 240 DF, p-value: 1.109e-15Anova(shanMod)## Anova Table (Type II tests)
##
## Response: Shannon
## Sum Sq Df F value Pr(>F)
## PopID 80.756 20 7.1663 1.319e-15 ***
## Species 10.113 2 8.9744 0.0001743 ***
## ReproductiveModeSimple 4.038 2 3.5829 0.0292918 *
## Residuals 135.227 240
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1summary(chaoMod)##
## Call:
## lm(formula = Chao1 ~ PopID + Species, data = alphaDiv)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.71014 -0.30002 0.02297 0.30643 1.05792
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.28383 0.23645 13.888 < 2e-16 ***
## PopIDM108 0.44896 0.23040 1.949 0.052496 .
## PopIDM109 0.61489 0.23581 2.608 0.009686 **
## PopIDM26 0.72523 0.24226 2.994 0.003042 **
## PopIDM28 0.48852 0.20582 2.373 0.018404 *
## PopIDM31 0.05747 0.25959 0.221 0.824972
## PopIDM34 -0.42703 0.27221 -1.569 0.118010
## PopIDM44 0.28525 0.21261 1.342 0.180971
## PopIDM67 -0.86959 0.23553 -3.692 0.000275 ***
## PopIDM70 -0.75293 0.26069 -2.888 0.004225 **
## PopIDM71 -0.83685 0.22182 -3.773 0.000203 ***
## PopIDM72 0.44506 0.20267 2.196 0.029045 *
## PopIDM78 -1.79411 0.40811 -4.396 1.65e-05 ***
## PopIDM79 0.18240 0.21528 0.847 0.397694
## PopIDM83 0.15393 0.21909 0.703 0.482985
## PopIDM84 0.42829 0.22182 1.931 0.054679 .
## PopIDM85 0.94649 0.21531 4.396 1.65e-05 ***
## PopIDM86 0.55613 0.38518 1.444 0.150081
## PopIDM89 0.36098 0.21234 1.700 0.090415 .
## PopIDM90 0.01105 0.19975 0.055 0.955916
## PopIDR12 -0.57231 0.23581 -2.427 0.015955 *
## SpeciesOLI 0.68218 0.17095 3.991 8.74e-05 ***
## SpeciesVUL 0.37232 0.17466 2.132 0.034044 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5267 on 242 degrees of freedom
## Multiple R-squared: 0.5293, Adjusted R-squared: 0.4865
## F-statistic: 12.37 on 22 and 242 DF, p-value: < 2.2e-16Anova(chaoMod)## Anova Table (Type II tests)
##
## Response: Chao1
## Sum Sq Df F value Pr(>F)
## PopID 67.399 20 12.1496 < 2.2e-16 ***
## Species 5.101 2 9.1957 0.0001415 ***
## Residuals 67.124 242
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1There is a clear effect of species, and population on alpha diversity in the samples. The populations M78 and M67 have the lowest alpha diversity, while M85 shows highest diversity in the models and plots. In the Shannon index, there is also some statistic difference for the reproductive mode. Sexually propagating animals show the highest alpha diversity. This coincides with the expression of periculin in the egg fleck and might be a consequence of the destroyed homoeostasis of the microbiota in these animals.
Finally we can have a look on the pairwise comparisons.
Anova(simpMod)## Anova Table (Type II tests)
##
## Response: Simpson
## Sum Sq Df F value Pr(>F)
## PopID 5.7527 20 6.0251 7.434e-13 ***
## Species 0.8197 2 8.5854 0.0002499 ***
## Residuals 11.5530 242
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova(shanMod)## Anova Table (Type II tests)
##
## Response: Shannon
## Sum Sq Df F value Pr(>F)
## PopID 80.756 20 7.1663 1.319e-15 ***
## Species 10.113 2 8.9744 0.0001743 ***
## ReproductiveModeSimple 4.038 2 3.5829 0.0292918 *
## Residuals 135.227 240
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova(chaoMod)## Anova Table (Type II tests)
##
## Response: Chao1
## Sum Sq Df F value Pr(>F)
## PopID 67.399 20 12.1496 < 2.2e-16 ***
## Species 5.101 2 9.1957 0.0001415 ***
## Residuals 67.124 242
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1TukeyHSD(aov(simpMod))## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = simpMod)
##
## $PopID
## diff lwr upr p adj
## M108-M107 0.019213394 -0.325400158 0.363826945 1.0000000
## M109-M107 0.268291064 -0.084432152 0.621014279 0.4220948
## M26-M107 0.180923214 -0.181465441 0.543311868 0.9671343
## M28-M107 0.210073212 -0.092960012 0.513106436 0.6074779
## M31-M107 0.172587803 -0.216094250 0.561269855 0.9909131
## M34-M107 0.192973661 -0.214316025 0.600263348 0.9810422
## M44-M107 0.188941744 -0.128998666 0.506882153 0.8478043
## M67-M107 -0.362616483 -0.715339698 -0.009893268 0.0361681
## M70-M107 -0.068445582 -0.442565048 0.305673884 1.0000000
## M71-M107 -0.099122430 -0.430872703 0.232627843 0.9999643
## M72-M107 0.034943603 -0.268089621 0.337976827 1.0000000
## M78-M107 -0.184197789 -0.795132319 0.426736740 0.9999588
## M79-M107 0.056710004 -0.265280765 0.378700773 1.0000000
## M83-M107 -0.033733892 -0.360292135 0.292824350 1.0000000
## M84-M107 0.205475656 -0.126274617 0.537225928 0.7960575
## M85-M107 0.260570771 -0.061419998 0.582561541 0.3029158
## M86-M107 0.289732328 -0.229462186 0.808926843 0.9075285
## M89-M107 0.042561710 -0.250921745 0.336045165 1.0000000
## M90-M107 0.121734394 -0.173375021 0.416843808 0.9963106
## R12-M107 -0.078146873 -0.430870088 0.274576342 0.9999998
## M109-M108 0.249077670 -0.095535882 0.593691221 0.5251314
## M26-M108 0.161709820 -0.192790349 0.516209990 0.9875813
## M28-M108 0.190859818 -0.102694211 0.484413847 0.7226320
## M31-M108 0.153374409 -0.227963461 0.534712279 0.9973342
## M34-M108 0.173760268 -0.226526782 0.574047317 0.9930626
## M44-M108 0.169728350 -0.139190632 0.478647332 0.9190472
## M67-M108 -0.381829876 -0.726443428 -0.037216325 0.0133830
## M70-M108 -0.087658976 -0.454142547 0.278824596 0.9999991
## M71-M108 -0.118335824 -0.441450473 0.204778826 0.9992486
## M72-M108 0.015730209 -0.277823821 0.309284238 1.0000000
## M78-M108 -0.203411183 -0.809699756 0.402877389 0.9997883
## M79-M108 0.037496610 -0.275589462 0.350582683 1.0000000
## M83-M108 -0.052947286 -0.370728845 0.264834272 1.0000000
## M84-M108 0.186262262 -0.136852388 0.509376912 0.8794462
## M85-M108 0.241357378 -0.071728695 0.554443450 0.3953962
## M86-M108 0.270518935 -0.243200617 0.784238486 0.9447680
## M89-M108 0.023348316 -0.260337023 0.307033655 1.0000000
## M90-M108 0.102521000 -0.182846132 0.387888132 0.9994273
## R12-M108 -0.097360267 -0.441973819 0.247253284 0.9999855
## M26-M109 -0.087367850 -0.449756504 0.275020805 0.9999990
## M28-M109 -0.058217852 -0.361251076 0.244815372 1.0000000
## M31-M109 -0.095703261 -0.484385313 0.292978792 0.9999985
## M34-M109 -0.075317402 -0.482607089 0.331972284 1.0000000
## M44-M109 -0.079349320 -0.397289729 0.238591090 0.9999981
## M67-M109 -0.630907546 -0.983630761 -0.278184331 0.0000001
## M70-M109 -0.336736645 -0.710856111 0.037382821 0.1409852
## M71-M109 -0.367413494 -0.699163766 -0.035663221 0.0134696
## M72-M109 -0.233347461 -0.536380685 0.069685763 0.3976013
## M78-M109 -0.452488853 -1.063423383 0.158445677 0.4757642
## M79-M109 -0.211581060 -0.533571829 0.110409710 0.7048733
## M83-M109 -0.302024956 -0.628583198 0.024533286 0.1109854
## M84-M109 -0.062815408 -0.394565681 0.268934865 1.0000000
## M85-M109 -0.007720292 -0.329711061 0.314270477 1.0000000
## M86-M109 0.021441265 -0.497753250 0.540635779 1.0000000
## M89-M109 -0.225729354 -0.519212809 0.067754101 0.3999205
## M90-M109 -0.146556670 -0.441666084 0.148552745 0.9688563
## R12-M109 -0.346437937 -0.699161152 0.006285278 0.0609831
## M28-M26 0.029149998 -0.285080817 0.343380813 1.0000000
## M31-M26 -0.008335411 -0.405809446 0.389138624 1.0000000
## M34-M26 0.012050448 -0.403637848 0.427738743 1.0000000
## M44-M26 0.008018530 -0.320611923 0.336648983 1.0000000
## M67-M26 -0.543539697 -0.905928351 -0.181151042 0.0000298
## M70-M26 -0.249368796 -0.632614468 0.133876876 0.7213473
## M71-M26 -0.280045644 -0.622054582 0.061963294 0.2818894
## M72-M26 -0.145979611 -0.460210426 0.168251204 0.9847417
## M78-M26 -0.365121003 -0.981686385 0.251444378 0.8517284
## M79-M26 -0.124213210 -0.456763846 0.208337427 0.9990118
## M83-M26 -0.214657106 -0.551632114 0.122317902 0.7551789
## M84-M26 0.024552442 -0.317456496 0.366561380 1.0000000
## M85-M26 0.079647558 -0.252903079 0.412198194 0.9999991
## M86-M26 0.108809115 -0.416999610 0.634617839 0.9999999
## M89-M26 -0.138361504 -0.443393320 0.166670312 0.9883579
## M90-M26 -0.059188820 -0.365785358 0.247407717 1.0000000
## R12-M26 -0.259070087 -0.621458742 0.103318567 0.5470988
## M31-M28 -0.037485409 -0.381708020 0.306737202 1.0000000
## M34-M28 -0.017099550 -0.382202764 0.348003664 1.0000000
## M44-M28 -0.021131468 -0.282859518 0.240596582 1.0000000
## M67-M28 -0.572689695 -0.875722918 -0.269656471 0.0000000
## M70-M28 -0.278518794 -0.606208957 0.049171370 0.2199190
## M71-M28 -0.309195642 -0.587536641 -0.030854643 0.0128919
## M72-M28 -0.175129609 -0.418531752 0.068272534 0.5342313
## M78-M28 -0.394271001 -0.977928836 0.189386833 0.6560621
## M79-M28 -0.153363208 -0.419996896 0.113270480 0.8815440
## M83-M28 -0.243807104 -0.515938973 0.028324764 0.1465833
## M84-M28 -0.004597556 -0.282938555 0.273743443 1.0000000
## M85-M28 0.050497560 -0.216136129 0.317131248 1.0000000
## M86-M28 0.079659117 -0.407145169 0.566463402 1.0000000
## M89-M28 -0.167511502 -0.398915924 0.063892920 0.5220242
## M90-M28 -0.088338818 -0.321801952 0.145124316 0.9988192
## R12-M28 -0.288220085 -0.591253309 0.014813139 0.0850220
## M34-M31 0.020385859 -0.418413595 0.459185312 1.0000000
## M44-M31 0.016353941 -0.341062029 0.373769911 1.0000000
## M67-M31 -0.535204286 -0.923886338 -0.146522233 0.0002452
## M70-M31 -0.241033385 -0.649231188 0.167164418 0.8548950
## M71-M31 -0.271710233 -0.641464631 0.098044165 0.4917075
## M72-M31 -0.137644200 -0.481866812 0.206578411 0.9975287
## M78-M31 -0.356785592 -0.989162909 0.275591725 0.8986697
## M79-M31 -0.115877799 -0.476901519 0.245145922 0.9998913
## M83-M31 -0.206321695 -0.571424909 0.158781519 0.8973185
## M84-M31 0.032887853 -0.336866545 0.402642251 1.0000000
## M85-M31 0.087982969 -0.273040752 0.449006689 0.9999988
## M86-M31 0.117144526 -0.427119212 0.661408263 0.9999999
## M89-M31 -0.130026093 -0.465872203 0.205820017 0.9983889
## M90-M31 -0.050853409 -0.388121312 0.286414493 1.0000000
## R12-M31 -0.250734676 -0.639416729 0.137947376 0.7353870
## M44-M34 -0.004031918 -0.381599562 0.373535727 1.0000000
## M67-M34 -0.555590144 -0.962879831 -0.148300458 0.0003028
## M70-M34 -0.261419243 -0.687372993 0.164534506 0.8083787
## M71-M34 -0.292096092 -0.681363950 0.097171767 0.4493933
## M72-M34 -0.158030059 -0.523133273 0.207073155 0.9933006
## M78-M34 -0.377171451 -1.021152989 0.266810087 0.8635955
## M79-M34 -0.136263657 -0.517248273 0.244720959 0.9994625
## M83-M34 -0.226707554 -0.611560133 0.158145025 0.8574778
## M84-M34 0.012501994 -0.376765864 0.401769853 1.0000000
## M85-M34 0.067597110 -0.313387506 0.448581726 1.0000000
## M86-M34 0.096758667 -0.460945705 0.654463039 1.0000000
## M89-M34 -0.150411952 -0.507628638 0.206804734 0.9952237
## M90-M34 -0.071239268 -0.429793014 0.287314478 1.0000000
## R12-M34 -0.271120535 -0.678410221 0.136169152 0.6824325
## M67-M44 -0.551558226 -0.869498636 -0.233617817 0.0000004
## M70-M44 -0.257387326 -0.598910110 0.084135459 0.4406173
## M71-M44 -0.288064174 -0.582564906 0.006436558 0.0637383
## M72-M44 -0.153998141 -0.415726192 0.107729909 0.8587294
## M78-M44 -0.373139533 -0.964674348 0.218395281 0.7699602
## M79-M44 -0.132231740 -0.415693365 0.151229885 0.9840219
## M83-M44 -0.222675636 -0.511315071 0.065963798 0.3939389
## M84-M44 0.016533912 -0.277966820 0.311034644 1.0000000
## M85-M44 0.071629028 -0.211832598 0.355090653 0.9999977
## M86-M44 0.100790585 -0.395430517 0.597011686 0.9999999
## M89-M44 -0.146380034 -0.396989244 0.104229175 0.8664106
## M90-M44 -0.067207350 -0.319718744 0.185304044 0.9999945
## R12-M44 -0.267088617 -0.585029027 0.050851792 0.2385305
## M70-M67 0.294170901 -0.079948565 0.668290367 0.3567695
## M71-M67 0.263494053 -0.068256220 0.595244325 0.3375773
## M72-M67 0.397560085 0.094526861 0.700593309 0.0007001
## M78-M67 0.178418693 -0.432515836 0.789353223 0.9999752
## M79-M67 0.419326487 0.097335718 0.741317256 0.0008146
## M83-M67 0.328882590 0.002324348 0.655440833 0.0461147
## M84-M67 0.568092138 0.236341866 0.899842411 0.0000006
## M85-M67 0.623187254 0.301196485 0.945178023 0.0000000
## M86-M67 0.652348811 0.133154296 1.171543326 0.0016478
## M89-M67 0.405178192 0.111694737 0.698661647 0.0002310
## M90-M67 0.484350876 0.189241462 0.779460291 0.0000022
## R12-M67 0.284469609 -0.068253606 0.637192825 0.3091257
## M71-M70 -0.030676848 -0.385091791 0.323738095 1.0000000
## M72-M70 0.103389184 -0.224300979 0.431079348 0.9999165
## M78-M70 -0.115752208 -0.739284651 0.507780236 1.0000000
## M79-M70 0.125155586 -0.220141046 0.470452217 0.9993506
## M83-M70 0.034711690 -0.314848025 0.384271404 1.0000000
## M84-M70 0.273921238 -0.080493705 0.628336181 0.3903185
## M85-M70 0.329016353 -0.016280278 0.674312985 0.0834989
## M86-M70 0.358177910 -0.175783373 0.892139194 0.6687808
## M89-M70 0.111007291 -0.207872383 0.429886966 0.9996358
## M90-M70 0.190179975 -0.130196795 0.510556746 0.8490435
## R12-M70 -0.009701291 -0.383820757 0.364418175 1.0000000
## M72-M71 0.134066033 -0.144274966 0.412407032 0.9773652
## M78-M71 -0.085075359 -0.684145953 0.513995234 1.0000000
## M79-M71 0.155832434 -0.143036488 0.454701356 0.9496754
## M83-M71 0.065388538 -0.238395699 0.369172775 0.9999999
## M84-M71 0.304598086 -0.004760639 0.613956810 0.0593734
## M85-M71 0.359693201 0.060824280 0.658562123 0.0035989
## M86-M71 0.388854758 -0.116325923 0.894035440 0.3983862
## M89-M71 0.141684140 -0.126228375 0.409596654 0.9425275
## M90-M71 0.220856824 -0.048835860 0.490549507 0.2816857
## R12-M71 0.020975557 -0.310774716 0.352725829 1.0000000
## M78-M72 -0.219141392 -0.802799227 0.364516442 0.9989384
## M79-M72 0.021766401 -0.244867287 0.288400090 1.0000000
## M83-M72 -0.068677495 -0.340809363 0.203454374 0.9999978
## M84-M72 0.170532053 -0.107808946 0.448873052 0.8106371
## M85-M72 0.225627169 -0.041006519 0.492260857 0.2268384
## M86-M72 0.254788726 -0.232015560 0.741593011 0.9478359
## M89-M72 0.007618107 -0.223786315 0.239022529 1.0000000
## M90-M72 0.086790791 -0.146672343 0.320253925 0.9990744
## R12-M72 -0.113090476 -0.416123700 0.189942748 0.9990235
## M79-M78 0.240907794 -0.352813849 0.834629436 0.9970124
## M83-M78 0.150463897 -0.445747155 0.746674949 0.9999978
## M84-M78 0.389673445 -0.209397149 0.988744039 0.7218870
## M85-M78 0.444768561 -0.148953081 1.038490203 0.4527683
## M86-M78 0.473930118 -0.246063130 1.193923366 0.7019185
## M89-M78 0.226759499 -0.351997679 0.805516677 0.9981107
## M90-M78 0.305932183 -0.273651201 0.885515568 0.9435248
## R12-M78 0.106050916 -0.504883613 0.716985446 1.0000000
## M83-M79 -0.090443896 -0.383538889 0.202651096 0.9999408
## M84-M79 0.148765652 -0.150103270 0.447634573 0.9681169
## M85-M79 0.203860767 -0.084136532 0.491858067 0.5667594
## M86-M79 0.233022324 -0.265803630 0.731848279 0.9837722
## M89-M79 -0.014148294 -0.269876524 0.241579935 1.0000000
## M90-M79 0.065024390 -0.192568226 0.322617005 0.9999978
## R12-M79 -0.134856877 -0.456847647 0.187133892 0.9955318
## M84-M83 0.239209548 -0.064574689 0.542993785 0.3540089
## M85-M83 0.294304664 0.001209671 0.587399656 0.0477141
## M86-M83 0.323466221 -0.178320151 0.825252593 0.7365308
## M89-M83 0.076295602 -0.185160240 0.337751444 0.9999755
## M90-M83 0.155468286 -0.107811386 0.418747958 0.8548499
## R12-M83 -0.044412981 -0.370971223 0.282145261 1.0000000
## M85-M84 0.055095116 -0.243773806 0.353964037 1.0000000
## M86-M84 0.084256673 -0.420924009 0.589437354 1.0000000
## M89-M84 -0.162913946 -0.430826460 0.104998568 0.8203072
## M90-M84 -0.083741262 -0.353433946 0.185951421 0.9999348
## R12-M84 -0.283622529 -0.615372802 0.048127744 0.2109243
## M86-M85 0.029161557 -0.469664398 0.527987512 1.0000000
## M89-M85 -0.218009062 -0.473737291 0.037719168 0.2152699
## M90-M85 -0.138836378 -0.396428993 0.118756238 0.9316658
## R12-M85 -0.338717645 -0.660708414 -0.016726876 0.0271592
## M89-M86 -0.247170619 -0.728088300 0.233747062 0.9561666
## M90-M86 -0.167997935 -0.649909591 0.313913721 0.9996283
## R12-M86 -0.367879202 -0.887073716 0.151315313 0.5648014
## M90-M89 0.079172684 -0.141753635 0.300099003 0.9994472
## R12-M89 -0.120708583 -0.414192038 0.172774872 0.9964458
## R12-M90 -0.199881267 -0.494990681 0.095228148 0.6511985
##
## $Species
## diff lwr upr p adj
## OLI-CIR 0.16630613 0.03973661 0.2928756 0.0061370
## VUL-CIR 0.20972749 0.06229557 0.3571594 0.0026458
## VUL-OLI 0.04342137 -0.04744624 0.1342890 0.4984479TukeyHSD(aov(shanMod))## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = shanMod)
##
## $PopID
## diff lwr upr p adj
## M108-M107 0.268667819 -0.915352007 1.45268764 0.9999996
## M109-M107 1.154821697 -0.057061231 2.36670462 0.0834317
## M26-M107 1.232022859 -0.013068489 2.47711421 0.0562408
## M28-M107 0.736880899 -0.304277649 1.77803945 0.5668822
## M31-M107 0.722758564 -0.612671381 2.05818851 0.9289949
## M34-M107 0.481000927 -0.918360942 1.88036280 0.9996965
## M44-M107 0.775958473 -0.316418035 1.86833498 0.5596522
## M67-M107 -1.047537142 -2.259420070 0.16434579 0.1944955
## M70-M107 -0.256454693 -1.541850648 1.02894126 1.0000000
## M71-M107 -0.417760280 -1.557584580 0.72206402 0.9992386
## M72-M107 0.527413943 -0.513744605 1.56857249 0.9619185
## M78-M107 -0.633481997 -2.732524801 1.46556081 0.9999581
## M79-M107 0.264682481 -0.841610214 1.37097517 0.9999991
## M83-M107 0.086633488 -1.035352085 1.20861906 1.0000000
## M84-M107 0.582482032 -0.557342268 1.72230633 0.9584811
## M85-M107 1.400606061 0.294313366 2.50689875 0.0014275
## M86-M107 0.672486082 -1.111357287 2.45632945 0.9988747
## M89-M107 0.219014597 -0.789332950 1.22736214 0.9999998
## M90-M107 0.410004549 -0.603929453 1.42393855 0.9971351
## R12-M107 -0.235372697 -1.447255625 0.97651023 1.0000000
## M109-M108 0.886153878 -0.297865947 2.07017370 0.4544826
## M26-M108 0.963355040 -0.254633130 2.18134321 0.3453938
## M28-M108 0.468213081 -0.540376945 1.47680311 0.9848373
## M31-M108 0.454090745 -0.856106131 1.76428762 0.9996575
## M34-M108 0.212333109 -1.162969170 1.58763539 1.0000000
## M44-M108 0.507290655 -0.554090122 1.56867143 0.9791278
## M67-M108 -1.316204961 -2.500224786 -0.13218514 0.0127654
## M70-M108 -0.525122512 -1.784283136 0.73403811 0.9957549
## M71-M108 -0.686428099 -1.796582210 0.42372601 0.7982246
## M72-M108 0.258746124 -0.749843901 1.26733615 0.9999970
## M78-M108 -0.902149816 -2.985230085 1.18093045 0.9932401
## M79-M108 -0.003985338 -1.079683364 1.07171269 1.0000000
## M83-M108 -0.182034331 -1.273865061 0.90979640 1.0000000
## M84-M108 0.313814214 -0.796339898 1.42396833 0.9999853
## M85-M108 1.131938242 0.056240216 2.20763627 0.0270541
## M86-M108 0.403818263 -1.361214281 2.16885081 0.9999996
## M89-M108 -0.049653222 -1.024336499 0.92503006 1.0000000
## M90-M108 0.141336731 -0.839124833 1.12179829 1.0000000
## R12-M108 -0.504040516 -1.688060341 0.67997931 0.9945188
## M26-M109 0.077201162 -1.167890186 1.32229251 1.0000000
## M28-M109 -0.417940797 -1.459099346 0.62321775 0.9973948
## M31-M109 -0.432063133 -1.767493078 0.90336681 0.9998770
## M34-M109 -0.673820770 -2.073182639 0.72554110 0.9773982
## M44-M109 -0.378863223 -1.471239732 0.71351329 0.9996540
## M67-M109 -2.202358839 -3.414241767 -0.99047591 0.0000001
## M70-M109 -1.411276390 -2.696672345 -0.12588044 0.0152573
## M71-M109 -1.572581977 -2.712406277 -0.43275768 0.0002351
## M72-M109 -0.627407754 -1.668566302 0.41375079 0.8313677
## M78-M109 -1.788303694 -3.887346498 0.31073911 0.2161813
## M79-M109 -0.890139216 -1.996431911 0.21615348 0.3133016
## M83-M109 -1.068188209 -2.190173782 0.05379736 0.0841867
## M84-M109 -0.572339665 -1.712163965 0.56748464 0.9651186
## M85-M109 0.245784364 -0.860508331 1.35207706 0.9999998
## M86-M109 -0.482335615 -2.266178984 1.30150775 0.9999929
## M89-M109 -0.935807100 -1.944154647 0.07254045 0.1075074
## M90-M109 -0.744817148 -1.758751150 0.26911685 0.4922593
## R12-M109 -1.390194394 -2.602077322 -0.17831147 0.0079417
## M28-M26 -0.495141960 -1.574773081 0.58448916 0.9867739
## M31-M26 -0.509264295 -1.874901645 0.85637306 0.9990311
## M34-M26 -0.751021932 -2.179239659 0.67719580 0.9454193
## M44-M26 -0.456064386 -1.585169639 0.67304087 0.9971764
## M67-M26 -2.279560001 -3.524651349 -1.03446865 0.0000001
## M70-M26 -1.488477553 -2.805229235 -0.17172587 0.0099540
## M71-M26 -1.649783139 -2.824854060 -0.47471222 0.0001568
## M72-M26 -0.704608916 -1.784240037 0.37502220 0.7161787
## M78-M26 -1.865504857 -3.983894086 0.25288437 0.1685547
## M79-M26 -0.967340379 -2.109914560 0.17523380 0.2259716
## M83-M26 -1.145389371 -2.303164765 0.01238602 0.0563672
## M84-M26 -0.649540827 -1.824611748 0.52553009 0.9146032
## M85-M26 0.168583201 -0.973990980 1.31115738 1.0000000
## M86-M26 -0.559536777 -2.366105182 1.24703163 0.9999371
## M89-M26 -1.013008262 -2.061033553 0.03501703 0.0722077
## M90-M26 -0.822018310 -1.875419656 0.23138304 0.3713590
## R12-M26 -1.467395556 -2.712486904 -0.22230421 0.0051405
## M31-M28 -0.014122335 -1.196798973 1.16855430 1.0000000
## M34-M28 -0.255879972 -1.510297978 0.99853803 0.9999999
## M44-M28 0.039077574 -0.860165068 0.93832022 1.0000000
## M67-M28 -1.784418041 -2.825576590 -0.74325949 0.0000005
## M70-M28 -0.993335593 -2.119210209 0.13253902 0.1660609
## M71-M28 -1.154641180 -2.110962426 -0.19831993 0.0034074
## M72-M28 -0.209466957 -1.045745627 0.62681171 0.9999980
## M78-M28 -1.370362897 -3.375688702 0.63496291 0.6345224
## M79-M28 -0.472198419 -1.388295803 0.44389897 0.9548575
## M83-M28 -0.650247412 -1.585235389 0.28474057 0.6011272
## M84-M28 -0.154398867 -1.110720114 0.80192238 1.0000000
## M85-M28 0.663725161 -0.252372223 1.57982255 0.5201264
## M86-M28 -0.064394818 -1.736952158 1.60816252 1.0000000
## M89-M28 -0.517866302 -1.312923324 0.27719072 0.7194437
## M90-M28 -0.326876350 -1.129006672 0.47525397 0.9968349
## R12-M28 -0.972253596 -2.013412145 0.06890495 0.1014822
## M34-M31 -0.241757637 -1.749380450 1.26586518 1.0000000
## M44-M31 0.053199909 -1.174806338 1.28120616 1.0000000
## M67-M31 -1.770295706 -3.105725651 -0.43486576 0.0005647
## M70-M31 -0.979213258 -2.381695224 0.42326871 0.5935126
## M71-M31 -1.140518844 -2.410917341 0.12987965 0.1440394
## M72-M31 -0.195344621 -1.378021259 0.98733202 1.0000000
## M78-M31 -1.356240562 -3.528956281 0.81647516 0.7850136
## M79-M31 -0.458076084 -1.698477806 0.78232564 0.9991537
## M83-M31 -0.636125076 -1.890543082 0.61829293 0.9615119
## M84-M31 -0.140276532 -1.410675029 1.13012197 1.0000000
## M85-M31 0.677847496 -0.562554226 1.91824922 0.9226687
## M86-M31 -0.050272482 -1.920248437 1.81970347 1.0000000
## M89-M31 -0.503743967 -1.657640705 0.65015277 0.9925576
## M90-M31 -0.312754015 -1.471535732 0.84602770 0.9999931
## R12-M31 -0.958131261 -2.293561206 0.37729868 0.5397444
## M44-M34 0.294957546 -1.002285626 1.59220072 0.9999996
## M67-M34 -1.528538069 -2.927899938 -0.12917620 0.0164009
## M70-M34 -0.737455621 -2.200943294 0.72603205 0.9638644
## M71-M34 -0.898761207 -2.236203859 0.43868144 0.6653679
## M72-M34 0.046413015 -1.208004990 1.30083102 1.0000000
## M78-M34 -1.114482925 -3.327068313 1.09810246 0.9640084
## M79-M34 -0.216318447 -1.525301615 1.09266472 1.0000000
## M83-M34 -0.394367440 -1.716640118 0.92790524 0.9999652
## M84-M34 0.101481105 -1.235961546 1.43892376 1.0000000
## M85-M34 0.919605133 -0.389378035 2.22858830 0.5814412
## M86-M34 0.191485155 -1.724670000 2.10764031 1.0000000
## M89-M34 -0.261986330 -1.489307879 0.96533522 0.9999999
## M90-M34 -0.070996378 -1.302911786 1.16091903 1.0000000
## R12-M34 -0.716373624 -2.115735493 0.68298824 0.9577603
## M67-M44 -1.823495615 -2.915872124 -0.73111911 0.0000013
## M70-M44 -1.032413167 -2.205813766 0.14098743 0.1697352
## M71-M44 -1.193718754 -2.205561450 -0.18187606 0.0050542
## M72-M44 -0.248544531 -1.147787173 0.65069811 0.9999898
## M78-M44 -1.409440471 -3.441829926 0.62294898 0.6065967
## M79-M44 -0.511275993 -1.485190636 0.46263865 0.9462644
## M83-M44 -0.689324986 -1.681029495 0.30237952 0.6021533
## M84-M44 -0.193476441 -1.205319137 0.81836625 1.0000000
## M85-M44 0.624647587 -0.349267056 1.59856223 0.7443888
## M86-M44 -0.103472392 -1.808383936 1.60143915 1.0000000
## M89-M44 -0.556943876 -1.417984515 0.30409676 0.7308797
## M90-M44 -0.365953924 -1.233530070 0.50162222 0.9951070
## R12-M44 -1.011331170 -2.103707679 0.08104534 0.1099401
## M70-M67 0.791082449 -0.494313506 2.07647840 0.8045262
## M71-M67 0.629776862 -0.510047438 1.76960116 0.9149291
## M72-M67 1.574951085 0.533792537 2.61610963 0.0000238
## M78-M67 0.414055145 -1.684987659 2.51309795 1.0000000
## M79-M67 1.312219623 0.205926928 2.41851232 0.0046177
## M83-M67 1.134170630 0.012185057 2.25615620 0.0441804
## M84-M67 1.630019174 0.490194874 2.76984347 0.0001009
## M85-M67 2.448143203 1.341850508 3.55443590 0.0000000
## M86-M67 1.720023224 -0.063820145 3.50386659 0.0740531
## M89-M67 1.266551739 0.258204192 2.27489929 0.0016601
## M90-M67 1.457541691 0.443607689 2.47147569 0.0000888
## R12-M67 0.812164445 -0.399718483 2.02404737 0.6702812
## M71-M70 -0.161305587 -1.379000936 1.05638976 1.0000000
## M72-M70 0.783868636 -0.342005980 1.90974325 0.5989881
## M78-M70 -0.377027304 -2.519353894 1.76529929 1.0000000
## M79-M70 0.521137174 -0.665229571 1.70750392 0.9919860
## M83-M70 0.343088181 -0.857925622 1.54410198 0.9999825
## M84-M70 0.838936726 -0.378758624 2.05663207 0.6191409
## M85-M70 1.657060754 0.470694009 2.84342750 0.0001770
## M86-M70 0.928940775 -0.905638112 2.76351966 0.9620775
## M89-M70 0.475469291 -0.620134335 1.57107292 0.9930695
## M90-M70 0.666459243 -0.434288090 1.76720658 0.8255764
## R12-M70 0.021081997 -1.264313958 1.30647795 1.0000000
## M72-M71 0.945174223 -0.011147024 1.90149547 0.0569698
## M78-M71 -0.215721717 -2.274002528 1.84255909 1.0000000
## M79-M71 0.682442761 -0.344408119 1.70929364 0.6850927
## M83-M71 0.504393768 -0.539345103 1.54813264 0.9765164
## M84-M71 1.000242312 -0.062649328 2.06313395 0.0942988
## M85-M71 1.818366341 0.791515461 2.84521722 0.0000002
## M86-M71 1.090246362 -0.645448418 2.82594114 0.7760160
## M89-M71 0.636774877 -0.283716285 1.55726604 0.6113393
## M90-M71 0.827764829 -0.098842620 1.75437228 0.1501077
## R12-M71 0.182387583 -0.957436717 1.32221188 1.0000000
## M78-M72 -1.160895940 -3.166221745 0.84442986 0.8752333
## M79-M72 -0.262731462 -1.178828846 0.65336592 0.9999814
## M83-M72 -0.440780455 -1.375768433 0.49420752 0.9820216
## M84-M72 0.055068090 -0.901253157 1.01138934 1.0000000
## M85-M72 0.873192118 -0.042905266 1.78928950 0.0832143
## M86-M72 0.145072139 -1.527485202 1.81762948 1.0000000
## M89-M72 -0.308399346 -1.103456367 0.48665768 0.9983433
## M90-M72 -0.117409393 -0.919539715 0.68472093 1.0000000
## R12-M72 -0.762786640 -1.803945188 0.27837191 0.4976327
## M79-M78 0.898164478 -1.141738458 2.93806741 0.9917606
## M83-M78 0.720115485 -1.328340540 2.76857151 0.9995792
## M84-M78 1.215964030 -0.842316781 3.27424484 0.8542333
## M85-M78 2.034088058 -0.005814878 4.07399099 0.0516329
## M86-M78 1.305968079 -1.167777588 3.77971375 0.9433674
## M89-M78 0.852496595 -1.135991584 2.84098477 0.9940237
## M90-M78 1.043486547 -0.947840303 3.03481340 0.9471748
## R12-M78 0.398109300 -1.700933503 2.49715210 1.0000000
## M83-M79 -0.178048993 -1.185061863 0.82896388 1.0000000
## M84-M79 0.317799552 -0.709051328 1.34465043 0.9999378
## M85-M79 1.135923580 0.146425313 2.12542185 0.0078509
## M86-M79 0.407803601 -1.306057671 2.12166487 0.9999992
## M89-M79 -0.045667883 -0.924296400 0.83296063 1.0000000
## M90-M79 0.145322069 -0.739712087 1.03035622 1.0000000
## R12-M79 -0.500055178 -1.606347872 0.60623752 0.9887970
## M84-M83 0.495848544 -0.547890326 1.53958742 0.9804490
## M85-M83 1.313972573 0.306959702 2.32098544 0.0007844
## M86-M83 0.585852594 -1.138180051 2.30988524 0.9997435
## M89-M83 0.132381109 -0.765926283 1.03068850 1.0000000
## M90-M83 0.323371062 -0.581202626 1.22794475 0.9994648
## R12-M83 -0.322006185 -1.443991758 0.79997939 0.9999812
## M85-M84 0.818124028 -0.208726852 1.84497491 0.3315983
## M86-M84 0.090004050 -1.645690731 1.82569883 1.0000000
## M89-M84 -0.363467435 -1.283958597 0.55702373 0.9979006
## M90-M84 -0.172477483 -1.099084933 0.75412997 1.0000000
## R12-M84 -0.817854729 -1.957679029 0.32196957 0.5395860
## M86-M85 -0.728119979 -2.441981251 0.98574129 0.9946538
## M89-M85 -1.181591463 -2.060219980 -0.30296295 0.0004153
## M90-M85 -0.990601511 -1.875635667 -0.10556736 0.0115438
## R12-M85 -1.635978758 -2.742271452 -0.52968606 0.0000434
## M89-M86 -0.453471485 -2.105803688 1.19886072 0.9999909
## M90-M86 -0.262481533 -1.918228825 1.39326576 1.0000000
## R12-M86 -0.907858779 -2.691702148 0.87598459 0.9601211
## M90-M89 0.190989952 -0.568066505 0.95004641 0.9999979
## R12-M89 -0.454387294 -1.462734841 0.55396025 0.9891869
## R12-M90 -0.645377246 -1.659311249 0.36855676 0.7562254
##
## $Species
## diff lwr upr p adj
## OLI-CIR 0.49329645 0.05844870 0.9281442 0.0216821
## VUL-CIR 0.49150464 -0.01501893 0.9980282 0.0593467
## VUL-OLI -0.00179181 -0.31398055 0.3103969 0.9998990
##
## $ReproductiveModeSimple
## diff lwr upr p adj
## NR-ASEX 0.07527977 -0.22845386 0.3790134 0.8285485
## SEX-ASEX 0.22784165 -0.02623834 0.4819216 0.0889673
## SEX-NR 0.15256187 -0.18618096 0.4913047 0.5385271TukeyHSD(aov(chaoMod))## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = chaoMod)
##
## $PopID
## diff lwr upr p adj
## M108-M107 0.479950567 -0.350709648 1.310610783 0.8772257
## M109-M107 0.645872725 -0.204335111 1.496080561 0.4246334
## M26-M107 0.756218971 -0.117286544 1.629724486 0.1923080
## M28-M107 0.386702299 -0.343732247 1.117136846 0.9419521
## M31-M107 0.044191388 -0.892692017 0.981074793 1.0000000
## M34-M107 -0.396043363 -1.377778809 0.585692083 0.9972282
## M44-M107 0.316236045 -0.450130942 1.082603032 0.9962950
## M67-M107 -0.869586527 -1.719794363 -0.019378691 0.0384702
## M70-M107 -0.954338319 -1.856119908 -0.052556729 0.0251345
## M71-M107 -0.805865845 -1.605520275 -0.006211415 0.0457733
## M72-M107 0.476048303 -0.254386243 1.206482849 0.7186610
## M78-M107 -1.763127646 -3.235730816 -0.290524477 0.0039359
## M79-M107 0.147247991 -0.628882026 0.923378008 1.0000000
## M83-M107 0.096387857 -0.690751647 0.883527361 1.0000000
## M84-M107 0.459273512 -0.340380918 1.258927942 0.8829182
## M85-M107 0.977481399 0.201351381 1.753611416 0.0015755
## M86-M107 -0.095063331 -1.346535380 1.156408718 1.0000000
## M89-M107 0.010412917 -0.697002763 0.717828597 1.0000000
## M90-M107 -0.129082098 -0.840417009 0.582252813 1.0000000
## R12-M107 -0.541320764 -1.391528600 0.308887073 0.7559491
## M109-M108 0.165922157 -0.664738058 0.996582373 1.0000000
## M26-M108 0.276268404 -0.578222622 1.130759430 0.9998786
## M28-M108 -0.093248268 -0.800834062 0.614337525 1.0000000
## M31-M108 -0.435759180 -1.354940089 0.483421729 0.9809226
## M34-M108 -0.875993930 -1.840850145 0.088862284 0.1309101
## M44-M108 -0.163714522 -0.908336163 0.580907118 0.9999998
## M67-M108 -1.349537095 -2.180197310 -0.518876879 0.0000030
## M70-M108 -1.434288886 -2.317664834 -0.550912939 0.0000031
## M71-M108 -1.285816413 -2.064655443 -0.506977382 0.0000018
## M72-M108 -0.003902265 -0.711488058 0.703683529 1.0000000
## M78-M108 -2.243078214 -3.704482717 -0.781673710 0.0000159
## M79-M108 -0.332702577 -1.087368617 0.421963464 0.9916513
## M83-M108 -0.383562711 -1.149546801 0.382421380 0.9661265
## M84-M108 -0.020677055 -0.799516086 0.758161976 1.0000000
## M85-M108 0.497530831 -0.257135209 1.252196872 0.6992385
## M86-M108 -0.575013899 -1.813289037 0.663261240 0.9848130
## M89-M108 -0.469537650 -1.153335847 0.214260546 0.6256382
## M90-M108 -0.609032665 -1.296884673 0.078819342 0.1614578
## R12-M108 -1.021271331 -1.851931546 -0.190611116 0.0024650
## M26-M109 0.110346247 -0.763159269 0.983851762 1.0000000
## M28-M109 -0.259170425 -0.989604972 0.471264121 0.9995205
## M31-M109 -0.601681337 -1.538564742 0.335202068 0.7424993
## M34-M109 -1.041916088 -2.023651534 -0.060180642 0.0242174
## M44-M109 -0.329636680 -1.096003667 0.436730307 0.9937919
## M67-M109 -1.515459252 -2.365667088 -0.665251416 0.0000001
## M70-M109 -1.600211044 -2.501992633 -0.698429454 0.0000002
## M71-M109 -1.451738570 -2.251393000 -0.652084140 0.0000001
## M72-M109 -0.169824422 -0.900258968 0.560610124 0.9999995
## M78-M109 -2.409000371 -3.881603540 -0.936397202 0.0000024
## M79-M109 -0.498624734 -1.274754751 0.277505284 0.7419226
## M83-M109 -0.549484868 -1.336624372 0.237654636 0.5940273
## M84-M109 -0.186599212 -0.986253642 0.613055217 0.9999994
## M85-M109 0.331608674 -0.444521343 1.107738691 0.9942796
## M86-M109 -0.740936056 -1.992408105 0.510535993 0.8519785
## M89-M109 -0.635459808 -1.342875488 0.071955873 0.1433780
## M90-M109 -0.774954823 -1.486289734 -0.063619912 0.0170089
## R12-M109 -1.187193488 -2.037401324 -0.336985652 0.0001777
## M28-M26 -0.369516672 -1.126942012 0.387908669 0.9740919
## M31-M26 -0.712027584 -1.670103277 0.246048110 0.4688492
## M34-M26 -1.152262334 -2.154241878 -0.150282790 0.0076334
## M44-M26 -0.439982926 -1.232117308 0.352151456 0.9111155
## M67-M26 -1.625805499 -2.499311014 -0.752299983 0.0000000
## M70-M26 -1.710557290 -2.634336785 -0.786777795 0.0000000
## M71-M26 -1.562084816 -2.386466844 -0.737702789 0.0000000
## M72-M26 -0.280170668 -1.037596009 0.477254672 0.9991365
## M78-M26 -2.519346618 -4.005522453 -1.033170782 0.0000008
## M79-M26 -0.608970980 -1.410554615 0.192612655 0.4245138
## M83-M26 -0.659831114 -1.472079302 0.152417073 0.2960066
## M84-M26 -0.296945459 -1.121327487 0.527436569 0.9994055
## M85-M26 0.221262427 -0.580321208 1.022846063 0.9999900
## M86-M26 -0.851282302 -2.118697313 0.416132709 0.6664356
## M89-M26 -0.745806054 -1.481058028 -0.010554080 0.0424354
## M90-M26 -0.885301069 -1.624324665 -0.146277473 0.0038995
## R12-M26 -1.297539735 -2.171045250 -0.424034220 0.0000385
## M31-M28 -0.342510912 -1.172228801 0.487206978 0.9962766
## M34-M28 -0.782745662 -1.662794381 0.097303057 0.1556335
## M44-M28 -0.070466254 -0.701338365 0.560405857 1.0000000
## M67-M28 -1.256288827 -1.986723373 -0.525854280 0.0000005
## M70-M28 -1.341040618 -2.130908517 -0.551172720 0.0000007
## M71-M28 -1.192568144 -1.863484286 -0.521652003 0.0000002
## M72-M28 0.089346004 -0.497353142 0.676045150 1.0000000
## M78-M28 -2.149829946 -3.556685075 -0.742974816 0.0000180
## M79-M28 -0.239454308 -0.882151022 0.403242405 0.9990456
## M83-M28 -0.290314443 -0.946264029 0.365635144 0.9912561
## M84-M28 0.072571213 -0.598344929 0.743487354 1.0000000
## M85-M28 0.590779099 -0.051917614 1.233475813 0.1173364
## M86-M28 -0.481765630 -1.655163922 0.691632662 0.9965246
## M89-M28 -0.376289382 -0.934069094 0.181490329 0.6584873
## M90-M28 -0.515784397 -1.078526449 0.046957654 0.1204262
## R12-M28 -0.928023063 -1.658457609 -0.197588517 0.0013314
## M34-M31 -0.440234751 -1.497921678 0.617452176 0.9958675
## M44-M31 0.272044657 -0.589474645 1.133563960 0.9999154
## M67-M31 -0.913777915 -1.850661320 0.023105490 0.0657153
## M70-M31 -0.998529707 -1.982454086 -0.014605327 0.0421960
## M71-M31 -0.850057233 -1.741317216 0.041202750 0.0827010
## M72-M31 0.431856915 -0.397860974 1.261574805 0.9505133
## M78-M31 -1.807319034 -3.331608129 -0.283029940 0.0046437
## M79-M31 0.103056603 -0.767158861 0.973272067 1.0000000
## M83-M31 0.052196469 -0.827852250 0.932245188 1.0000000
## M84-M31 0.415082125 -0.476177858 1.306342107 0.9843104
## M85-M31 0.933290011 0.063074547 1.803505475 0.0210765
## M86-M31 -0.139254719 -1.451153892 1.172644454 1.0000000
## M89-M31 -0.033778471 -0.843305551 0.775748610 1.0000000
## M90-M31 -0.173273486 -0.986227669 0.639680698 0.9999999
## R12-M31 -0.585512151 -1.522395556 0.351371254 0.7834814
## M44-M34 0.712279408 -0.197813708 1.622372524 0.3657738
## M67-M34 -0.473543164 -1.455278610 0.508192282 0.9770166
## M70-M34 -0.558294956 -1.585018461 0.468428550 0.9260670
## M71-M34 -0.409822482 -1.348117920 0.528472956 0.9925289
## M72-M34 0.872091666 -0.007957053 1.752140385 0.0553397
## M78-M34 -1.367084283 -2.919344318 0.185175751 0.1684780
## M79-M34 0.543291354 -0.375038067 1.461620775 0.8528054
## M83-M34 0.492431220 -0.435221582 1.420084021 0.9404971
## M84-M34 0.855316875 -0.082978563 1.793612313 0.1263885
## M85-M34 1.373524762 0.455195341 2.291854183 0.0000321
## M86-M34 0.300980032 -1.043316592 1.645276655 0.9999997
## M89-M34 0.406456280 -0.454582666 1.267495226 0.9817910
## M90-M34 0.266961265 -0.597300545 1.131223075 0.9999398
## R12-M34 -0.145277401 -1.127012847 0.836458046 1.0000000
## M67-M44 -1.185822572 -1.952189559 -0.419455585 0.0000127
## M70-M44 -1.270574364 -2.093784562 -0.447364166 0.0000136
## M71-M44 -1.122101890 -1.831969625 -0.412234155 0.0000069
## M72-M44 0.159812258 -0.471059853 0.790684369 0.9999976
## M78-M44 -2.079363692 -3.505205579 -0.653521804 0.0000618
## M79-M44 -0.168988054 -0.852247008 0.514270900 0.9999984
## M83-M44 -0.219848188 -0.915587789 0.475891413 0.9999145
## M84-M44 0.143037467 -0.566830268 0.852905202 1.0000000
## M85-M44 0.661245353 -0.022013600 1.344504307 0.0712976
## M86-M44 -0.411299376 -1.607396063 0.784797311 0.9996955
## M89-M44 -0.305823128 -0.909894265 0.298248009 0.9621934
## M90-M44 -0.445318143 -1.053974326 0.163338040 0.5005306
## R12-M44 -0.857556809 -1.623923796 -0.091189822 0.0115822
## M70-M67 -0.084751792 -0.986533381 0.817029798 1.0000000
## M71-M67 0.063720682 -0.735933748 0.863375112 1.0000000
## M72-M67 1.345634830 0.615200284 2.076069376 0.0000000
## M78-M67 -0.893541119 -2.366144288 0.579062050 0.8230199
## M79-M67 1.016834518 0.240704501 1.792964536 0.0007203
## M83-M67 0.965974384 0.178834880 1.753113888 0.0025490
## M84-M67 1.328860040 0.529205610 2.128514469 0.0000015
## M85-M67 1.847067926 1.070937909 2.623197943 0.0000000
## M86-M67 0.774523196 -0.476948853 2.025995245 0.7971179
## M89-M67 0.879999444 0.172583764 1.587415125 0.0019879
## M90-M67 0.740504429 0.029169519 1.451839340 0.0310015
## R12-M67 0.328265764 -0.521942072 1.178473600 0.9984472
## M71-M70 0.148472474 -0.705813121 1.002758069 1.0000000
## M72-M70 1.430386622 0.640518723 2.220254520 0.0000001
## M78-M70 -0.808789328 -2.311758643 0.694179988 0.9326269
## M79-M70 1.101586310 0.269279591 1.933893029 0.0005827
## M83-M70 1.050726176 0.208143676 1.893308675 0.0018988
## M84-M70 1.413611831 0.559326236 2.267897426 0.0000017
## M85-M70 1.931819717 1.099512999 2.764126436 0.0000000
## M86-M70 0.859274988 -0.427791040 2.146341016 0.6772477
## M89-M70 0.964751236 0.196120235 1.733382237 0.0016806
## M90-M70 0.825256221 0.053016604 1.597495838 0.0221082
## R12-M70 0.413017555 -0.488764034 1.314799145 0.9869993
## M72-M71 1.281914148 0.610998007 1.952830289 0.0000000
## M78-M71 -0.957261801 -2.401268012 0.486744409 0.6896892
## M79-M71 0.953113836 0.232716969 1.673510703 0.0005873
## M83-M71 0.902253702 0.170008906 1.634498498 0.0023685
## M84-M71 1.265139357 0.519457757 2.010820958 0.0000008
## M85-M71 1.783347244 1.062950376 2.503744111 0.0000000
## M86-M71 0.710802514 -0.506890441 1.928495469 0.8670571
## M89-M71 0.816278762 0.170499553 1.462057972 0.0014696
## M90-M71 0.676783747 0.026713599 1.326853895 0.0309693
## R12-M71 0.264545082 -0.535109348 1.064199511 0.9998283
## M78-M72 -2.239175949 -3.646031079 -0.832320820 0.0000056
## M79-M72 -0.328800312 -0.971497025 0.313896401 0.9580889
## M83-M72 -0.379660446 -1.035610033 0.276289140 0.8755448
## M84-M72 -0.016774791 -0.687690932 0.654141351 1.0000000
## M85-M72 0.501433096 -0.141263618 1.144129809 0.3718473
## M86-M72 -0.571111634 -1.744509926 0.602286658 0.9747216
## M89-M72 -0.465635386 -1.023415097 0.092144326 0.2488381
## M90-M72 -0.605130401 -1.167872453 -0.042388349 0.0203386
## R12-M72 -1.017369066 -1.747803613 -0.286934520 0.0001885
## M79-M78 1.910375637 0.479262597 3.341488677 0.0004861
## M83-M78 1.859515503 0.422401963 3.296629043 0.0009282
## M84-M78 2.222401159 0.778394948 3.666407369 0.0000148
## M85-M78 2.740609045 1.309496005 4.171722085 0.0000000
## M86-M78 1.668064315 -0.067415163 3.403543794 0.0765155
## M89-M78 1.773540564 0.378498029 3.168583098 0.0013151
## M90-M78 1.634045548 0.237011517 3.031079580 0.0058174
## R12-M78 1.221806883 -0.250796286 2.694410052 0.2591282
## M83-M79 -0.050860134 -0.757339460 0.655619191 1.0000000
## M84-M79 0.312025521 -0.408371346 1.032422389 0.9932454
## M85-M79 0.830233408 0.136041616 1.524425199 0.0040107
## M86-M79 -0.242311322 -1.444686775 0.960064131 1.0000000
## M89-M79 -0.136835074 -0.753245152 0.479575004 0.9999998
## M90-M79 -0.276330089 -0.897234103 0.344573925 0.9906668
## R12-M79 -0.688568754 -1.464698772 0.087561263 0.1588752
## M84-M83 0.362885655 -0.369359140 1.095130451 0.9695160
## M85-M83 0.881093542 0.174614216 1.587572867 0.0018951
## M86-M83 -0.191451188 -1.400962463 1.018060087 1.0000000
## M89-M83 -0.085974940 -0.716190917 0.544241038 1.0000000
## M90-M83 -0.225469955 -0.860082111 0.409142201 0.9995113
## R12-M83 -0.637708620 -1.424848124 0.149430883 0.3008788
## M85-M84 0.518207886 -0.202188981 1.238604754 0.5347069
## M86-M84 -0.554336843 -1.772029798 0.663356112 0.9878667
## M89-M84 -0.448860595 -1.094639805 0.196918614 0.6023924
## M90-M84 -0.588355610 -1.238425759 0.061714538 0.1344962
## R12-M84 -1.000594276 -1.800248706 -0.200939846 0.0017815
## M86-M85 -1.072544730 -2.274920182 0.129830723 0.1519716
## M89-M85 -0.967068482 -1.583478559 -0.350658404 0.0000086
## M90-M85 -1.106563497 -1.727467511 -0.485659483 0.0000001
## R12-M85 -1.518802162 -2.294932180 -0.742672145 0.0000000
## M89-M86 0.105476248 -1.053732909 1.264685405 1.0000000
## M90-M86 -0.034018767 -1.195623812 1.127586278 1.0000000
## R12-M86 -0.446257432 -1.697729482 0.805214617 0.9994851
## M90-M89 -0.139495015 -0.672018193 0.393028163 0.9999958
## R12-M89 -0.551733681 -1.259149361 0.155682000 0.3725295
## R12-M90 -0.412238666 -1.123573576 0.299096245 0.8743051
##
## $Species
## diff lwr upr p adj
## OLI-CIR 0.3610004 0.05591592 0.66608496 0.0156291
## VUL-CIR 0.2098717 -0.14549980 0.56524321 0.3463080
## VUL-OLI -0.1511287 -0.37015701 0.06789953 0.2361629# do some nice plots
models <- list("Simpson" = simpMod ,"Shannon" = shanMod,"Chao1" = chaoMod)
effecSizePlots <- list()
for(vrb in c("Species","PopID")){
for(n in names(models)){
mod <- models[[n]]
prd <- make_predictions(mod, pred = vrb, data = alphaDiv)
effecSizePlots[[paste(n,vrb,sep="_")]] <-
ggplot(prd, aes_string(y = vrb, x=n ,xmin = "ymin", xmax = "ymax")) +
geom_point() +
geom_errorbar() +
ggtitle(paste0("Effect size of ",n))
}
}
pp <- plot_grid(plotlist=effecSizePlots, ncol = 3, rel_heights =c(0.3,0.7))
print(pp)
This last part is has to be taken carefully, as the design is unbalanced and Anova calculates Type II anovas, so only main effects are considered, which is fine as there are no interactions.